Analysis_Production_Exp7
Andrea Hofmann
30.09.2021 updated 01 August, 2022
1 Setup
2 Import
data <- read.csv2("data/d.csv", header = TRUE)3 Analysis
Just a short information on the amount of pause values to be found in the ungrouped condition to find out, if it is reasonable to model pause data for this condition.
# how many not zero pauses do we have in the ungrouped condition ---------------------------------
grou <- subset(data, data$condition=="grouped")
colSums(grou != 0)## X.1 X Subj group exp rep
## 658 658 658 658 658 658
## trial condition f0_range1 f0_range2 s4reln1 s8reln2
## NA 658 NA NA 658 658
## pause1 pause3 meanf0 list L1 L1t
## 80 522 658 658 NA NA
## H1 H1t f0n1 slopen1 L2 L2t
## NA NA NA NA NA NA
## H2 H2t f0n2 slopen2 s4dur s8dur
## NA NA NA NA 658 658
## p1dur p2dur p3dur p4dur fp0dur fp1dur
## 80 53 522 141 1 0
## fp2dur fp3dur fp4dur c1dur c2dur utt_dur.ms.
## 7 0 3 658 658 658
## filename rating_match item pause2 pause4 accuracy
## 658 NA 658 53 141 658
## rate
## 428rm(grou)
un_grou <- subset(data, data$condition=="ungrouped")
colSums(un_grou != 0)## X.1 X Subj group exp rep
## 643 643 643 643 643 643
## trial condition f0_range1 f0_range2 s4reln1 s8reln2
## NA 643 NA NA 643 643
## pause1 pause3 meanf0 list L1 L1t
## 118 186 643 643 NA NA
## H1 H1t f0n1 slopen1 L2 L2t
## NA NA NA NA NA NA
## H2 H2t f0n2 slopen2 s4dur s8dur
## NA NA NA NA 643 643
## p1dur p2dur p3dur p4dur fp0dur fp1dur
## 118 38 186 65 1 0
## fp2dur fp3dur fp4dur c1dur c2dur utt_dur.ms.
## 6 1 4 643 643 643
## filename rating_match item pause2 pause4 accuracy
## 643 NA 643 38 65 643
## rate
## 509rm(un_grou)- non-zero pause values in ungrouped condition:
- data points overall: 643
- pause1: 118
- pause2: 38
- pause3: 186 (29% of data points overall in ungrouped
condition)
- pause4: 65
- pause1: 118
- data points overall: 643
- non-zero pause values in grouped condition:
- data points overall: 658
- pause1: 80
- pause2: 53
- pause3: 522 (79% of data points overall in grouped condition)
- pause4: 141
- pause1: 80
- data points overall: 658
I would suggest, it makes sense to include pause3 in both conditions.
3.1 Contrast coding
- assign scaled sum contrast to all factors - all have two levels
(0.5, -0.5)
- intercept estimates average effects
- slopes estimate difference between factor levels
# Contrast coding ---------------------------------
## we have 4 predictors / factors with two levels (condition, group, exp, accuracy)
## therefore it is best to use sum contrast with 0.5 / -0.5
## within each predictor this compares averages of each level against each other - mean of these contrasts is 0
## interactions between predictors estimate average effects (estimates are at mean which is 0 = averaged)
# Prepare data ---------------------------------
## Subset for patient data - only RHDP / LHDP =================================
## model stimuli (stimuli that had to be repeated) will be "removed"
pDat <- subset(data, data$group!="model_stimuli")
### Factorize predictors variables #############################
########################################### condition
class(pDat$condition) #character## [1] "character"pDat$condition <- as.factor(pDat$condition)
## check level order
levels(pDat$condition) # "grouped" "ungrouped"## [1] "grouped" "ungrouped"########################################### group
class(pDat$group) #character ## [1] "character"pDat$group <- as.factor(pDat$group)
## check level order
levels(pDat$group) # "LHDP", "RHDP" ## [1] "LHDP" "RHDP"########################################### exp
class(pDat$exp) # character## [1] "character"pDat$exp <- as.factor(pDat$exp)
## check level order
levels(pDat$exp) # "reading_aloud" "repetition" ## [1] "reading_aloud" "repetition"########################################### accuracy
class(pDat$accuracy) # "character"## [1] "character"## factorize "accuracy"
pDat$accuracy <- as.factor(as.character(pDat$accuracy))
## check level order
levels(pDat$accuracy) # "correct" "incorrect"## [1] "correct" "incorrect"### contrast code predictor variables #############################
########################################### condition
contrasts(pDat$condition) <- contr.sum(2)/2
contrasts(pDat$condition)## [,1]
## grouped 0.5
## ungrouped -0.5# [,1]
# grouped 0.5
# ungrouped -0.5
########################################### group
## define sum contrast for group
contrasts(pDat$group) <- contr.sum(2)/2
contrasts(pDat$group)## [,1]
## LHDP 0.5
## RHDP -0.5# [,1]
# LHDP 0.5
# RHDP -0.5
########################################### exp
## assign sum contrast to experiment
contrasts(pDat$exp) <- contr.sum(2)/2
contrasts(pDat$exp)## [,1]
## reading_aloud 0.5
## repetition -0.5# [,1]
# reading_aloud 0.5
# repetition -0.5
########################################### accuracy
## define sum contrast for accuracy
contrasts(pDat$accuracy) <- contr.sum(2)/2
contrasts(pDat$accuracy)## [,1]
## correct 0.5
## incorrect -0.5# [,1]
# correct 0.5
# incorrect -0.53.2 Find correct / interpretable / parsimonious model to analyze data
3.2.1 Maximally complex model based on our hyoptheses
- the random effects structure is reduced iteratively in the next step to avoid overparameterization as best as possible
# Max models based on hypotheses ---------------------------------
## FL name2 =================================
summary(FL2_h <- lmer(s8reln2 ~ condition*accuracy+condition*group*exp+
group*exp*accuracy+
(1+condition*accuracy+exp*accuracy+condition*exp|Subj)+
(1+group*accuracy|item),
REML=TRUE,pDat))## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: s8reln2 ~ condition * accuracy + condition * group * exp + group *
## exp * accuracy + (1 + condition * accuracy + exp * accuracy +
## condition * exp | Subj) + (1 + group * accuracy | item)
## Data: pDat
##
## REML criterion at convergence: 8222.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2221 -0.6159 -0.0027 0.5753 4.0519
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 18.894 4.347
## condition1 7.272 2.697 -0.72
## accuracy1 2.962 1.721 -0.47 0.29
## exp1 10.784 3.284 0.05 -0.07 0.47
## condition1:accuracy1 51.198 7.155 -0.14 0.18 -0.80 -0.43
## accuracy1:exp1 13.125 3.623 -0.02 0.19 -0.03 -0.10 0.05
## condition1:exp1 12.211 3.494 -0.13 -0.19 -0.30 -0.24 0.39
## item (Intercept) 7.809 2.794
## group1 5.637 2.374 -0.99
## accuracy1 3.373 1.837 0.83 -0.88
## group1:accuracy1 10.218 3.197 0.92 -0.96 0.98
## Residual 27.753 5.268
##
##
##
##
##
##
##
## -0.84
##
##
##
##
##
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 43.4670 0.9944 45.3883 43.711 < 2e-16 ***
## condition1 4.7019 1.2593 31.3445 3.734 0.000752 ***
## accuracy1 -0.1099 0.6756 28.6661 -0.163 0.871936
## group1 3.9785 1.6936 33.0783 2.349 0.024943 *
## exp1 1.7234 1.3869 31.5163 1.243 0.223160
## condition1:accuracy1 9.5060 1.6397 26.0327 5.797 4.14e-06 ***
## condition1:group1 -0.4442 1.5411 38.2497 -0.288 0.774703
## condition1:exp1 0.2030 2.2923 29.2689 0.089 0.930022
## group1:exp1 -4.2497 1.7859 35.0612 -2.380 0.022904 *
## accuracy1:group1 -1.7295 1.2535 35.4761 -1.380 0.176306
## accuracy1:exp1 -0.7491 1.3584 26.8565 -0.551 0.585889
## condition1:group1:exp1 -2.2148 2.1714 31.6666 -1.020 0.315482
## accuracy1:group1:exp1 4.3390 2.6116 29.5572 1.661 0.107204
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# Model failed to converge with max|grad| = 0.0250477 (tol = 0.002, component 1)
# reduction needed
## F0 name2 =================================
summary(F02_h <- lmer(f0_range2 ~ condition*accuracy+condition*group*exp+
group*exp*accuracy+
(1+condition*accuracy+exp*accuracy+condition*exp|Subj)+
(1+group*accuracy|item),
REML=TRUE,pDat))## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: f0_range2 ~ condition * accuracy + condition * group * exp +
## group * exp * accuracy + (1 + condition * accuracy + exp *
## accuracy + condition * exp | Subj) + (1 + group * accuracy | item)
## Data: pDat
##
## REML criterion at convergence: 5846.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4215 -0.4556 -0.0352 0.3921 6.9327
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 3.32418 1.8232
## condition1 1.75620 1.3252 0.38
## accuracy1 0.22160 0.4707 0.44 1.00
## exp1 1.84342 1.3577 -0.14 -0.02 -0.03
## condition1:accuracy1 6.70690 2.5898 0.54 -0.06 -0.03 -0.38
## accuracy1:exp1 1.25981 1.1224 0.72 0.64 0.67 -0.14 0.17
## condition1:exp1 2.33252 1.5273 -0.03 0.56 0.54 0.77 -0.57
## item (Intercept) 0.04047 0.2012
## group1 0.62260 0.7890 -0.82
## accuracy1 0.16225 0.4028 0.32 -0.80
## group1:accuracy1 0.45821 0.6769 0.95 -0.59 0.00
## Residual 5.42452 2.3291
##
##
##
##
##
##
##
## 0.24
##
##
##
##
##
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.2720 0.3446 29.9811 15.298 1.05e-15 ***
## condition1 2.0080 0.3113 27.0110 6.451 6.50e-07 ***
## accuracy1 0.3662 0.2373 52.0728 1.543 0.128914
## group1 0.9976 0.6488 36.0793 1.538 0.132854
## exp1 2.2673 0.3369 27.3297 6.730 2.98e-07 ***
## condition1:accuracy1 3.0841 0.6293 27.2214 4.900 3.90e-05 ***
## condition1:group1 0.1046 0.6575 31.1750 0.159 0.874670
## condition1:exp1 0.8512 0.4699 34.8730 1.811 0.078710 .
## group1:exp1 2.8110 0.7055 34.3688 3.984 0.000334 ***
## accuracy1:group1 -0.2812 0.4659 44.2023 -0.604 0.549238
## accuracy1:exp1 0.3828 0.4848 25.3706 0.790 0.437107
## condition1:group1:exp1 1.8709 0.9410 36.0065 1.988 0.054433 .
## accuracy1:group1:exp1 1.4779 0.9517 21.8607 1.553 0.134785
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
# reduction needed
## Pause3 =================================
summary(P_h <- lmer(pause3 ~ condition*accuracy+condition*group*exp+
group*exp*accuracy+
(1+condition*accuracy+exp*accuracy+condition*exp|Subj)+
(1+group*accuracy|item),REML=TRUE,pDat))## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: pause3 ~ condition * accuracy + condition * group * exp + group *
## exp * accuracy + (1 + condition * accuracy + exp * accuracy +
## condition * exp | Subj) + (1 + group * accuracy | item)
## Data: pDat
##
## REML criterion at convergence: 7906.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4980 -0.3586 -0.0742 0.3440 6.8107
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 7.89696 2.8102
## condition1 9.62927 3.1031 0.24
## accuracy1 5.12554 2.2640 -0.20 0.43
## exp1 14.10797 3.7561 0.27 0.38 -0.19
## condition1:accuracy1 54.71130 7.3967 0.78 0.07 0.08 -0.31
## accuracy1:exp1 11.64825 3.4130 0.33 0.33 0.70 -0.47 0.74
## condition1:exp1 31.55293 5.6172 0.32 0.74 0.00 0.72 0.03
## item (Intercept) 0.00000 0.0000
## group1 0.34719 0.5892 NaN
## accuracy1 0.06886 0.2624 NaN 1.00
## group1:accuracy1 6.89504 2.6258 NaN -1.00 -1.00
## Residual 22.89651 4.7850
##
##
##
##
##
##
##
## -0.02
##
##
##
##
##
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.9994 0.5498 28.5552 16.370 4.78e-16 ***
## condition1 8.4077 0.6865 29.6301 12.246 4.03e-13 ***
## accuracy1 2.0065 0.6163 27.1148 3.256 0.00303 **
## group1 2.0490 0.9164 31.5150 2.236 0.03257 *
## exp1 4.6392 0.8287 24.9723 5.598 8.03e-06 ***
## condition1:accuracy1 17.5473 1.5823 28.6336 11.090 7.05e-12 ***
## condition1:group1 1.1401 1.3838 31.1602 0.824 0.41628
## condition1:exp1 7.1036 1.2255 25.4400 5.797 4.52e-06 ***
## group1:exp1 -0.7602 1.5975 28.2742 -0.476 0.63781
## accuracy1:group1 -0.2328 1.3381 28.6351 -0.174 0.86309
## accuracy1:exp1 0.3461 1.1330 31.3823 0.306 0.76199
## condition1:group1:exp1 -0.8169 2.4463 26.0634 -0.334 0.74109
## accuracy1:group1:exp1 2.0963 2.3228 32.2679 0.902 0.37348
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
# reduction needed4 rePCA
4.1 Indicator variables for all parameters
Create model.matrix dependent on model specification to merge these cols for all indicator variables with data.frame to model data
# generate model matrix ---------------------------------
pDat_mm <- pDat
mm <- model.matrix(~condition*accuracy+condition*group*exp+group*exp*accuracy, pDat)
colnames(mm)## [1] "(Intercept)" "condition1" "accuracy1"
## [4] "group1" "exp1" "condition1:accuracy1"
## [7] "condition1:group1" "condition1:exp1" "group1:exp1"
## [10] "accuracy1:group1" "accuracy1:exp1" "condition1:group1:exp1"
## [13] "accuracy1:group1:exp1"## create new col in pDat from mm cols =================================
# I renamed the col to reflect the level comparisons
## condition1 - grpd.ungr
## accuracy1 - cor.inc
## group1 - lhd.rhd
## exp1 - rd.rpt
pDat_mm$grpd.ungr <-mm[,2]
pDat_mm$cor.inc <-mm[,3]
pDat_mm$lhd.rhd <-mm[,4]
pDat_mm$rd.rpt <-mm[,5]
pDat_mm$grpd.ungr_cor.inc <-mm[,6]
pDat_mm$grpd.ungr_lhd.rhd <-mm[,7]
pDat_mm$grpd.ungr_rd.rpt <-mm[,8]
pDat_mm$lhd.rhd_rd.rpt <-mm[,9]
pDat_mm$cor.inc_lhd.rhd <-mm[,10]
pDat_mm$cor.inc_rd.rpt <-mm[,11]
pDat_mm$grpd.ungr_lhd.rhd_rd.rpt <-mm[,12]
pDat_mm$cor.inc_lhd.rhd_rd.rpt <-mm[,13]
# Clean workspace ---------------------------------
rm(mm)4.2 Model reduction - final lengthening (FL2)
# all "REML = FALSE" to use ANOVA model comparison later on
# INITIAL model with indicator variables ---------------------------------
summary(m0_FL2 <- lmer(s8reln2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd+cor.inc+
cor.inc_lhd.rhd|item),
REML = FALSE, pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## s8reln2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + cor.inc_rd.rpt + grpd.ungr_rd.rpt | Subj) +
## (1 + lhd.rhd + cor.inc + cor.inc_lhd.rhd | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8358.2 8626.6 -4127.1 8254.2 1237
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2136 -0.6172 -0.0114 0.5835 4.0524
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 17.629 4.199
## cor.inc 2.609 1.615 -0.45
## grpd.ungr 6.689 2.586 -0.73 0.26
## rd.rpt 9.973 3.158 0.05 0.49 -0.06
## grpd.ungr_cor.inc 49.570 7.041 -0.14 -0.81 0.19 -0.44
## cor.inc_rd.rpt 10.508 3.242 0.01 -0.03 0.21 -0.06 0.03
## grpd.ungr_rd.rpt 10.963 3.311 -0.14 -0.30 -0.20 -0.27 0.40
## item (Intercept) 6.752 2.598
## lhd.rhd 4.961 2.227 -1.00
## cor.inc 2.836 1.684 0.88 -0.92
## cor.inc_lhd.rhd 8.860 2.977 0.95 -0.97 0.99
## Residual 27.763 5.269
##
##
##
##
##
##
##
## -0.86
##
##
##
##
##
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 43.4630 0.9498 48.5381 45.761 < 2e-16 ***
## grpd.ungr 4.6890 1.1872 35.8285 3.950 0.000351 ***
## lhd.rhd 3.9835 1.6340 35.0837 2.438 0.019978 *
## cor.inc -0.1039 0.6432 30.5962 -0.162 0.872692
## rd.rpt 1.7185 1.3072 35.9554 1.315 0.196956
## grpd.ungr_lhd.rhd -0.4688 1.4848 40.3254 -0.316 0.753815
## grpd.ungr_cor.inc 9.5546 1.6008 26.7443 5.969 2.39e-06 ***
## cor.inc_lhd.rhd -1.7708 1.1995 39.0258 -1.476 0.147889
## grpd.ungr_rd.rpt 0.2958 2.1587 34.1389 0.137 0.891830
## lhd.rhd_rd.rpt -4.2403 1.7142 37.2433 -2.474 0.018056 *
## cor.inc_rd.rpt -0.7354 1.2808 29.9248 -0.574 0.570159
## grpd.ungr_lhd.rhd_rd.rpt -2.1463 2.0879 34.7450 -1.028 0.311067
## cor.inc_lhd.rhd_rd.rpt 4.3161 2.4791 34.0037 1.741 0.090731 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
# correlations seem okay for Subj
# but not perfect for items
'item (Intercept) 6.750 2.598
lhd.rhd 4.960 2.227 -1.00 #not optimal
cor.inc 2.834 1.683 0.88 -0.92
cor.inc_lhd.rhd 8.860 2.977 0.95 -0.97 0.99 #not optimal'## [1] "item (Intercept) 6.750 2.598 \n lhd.rhd 4.960 2.227 -1.00 #not optimal \n cor.inc 2.834 1.683 0.88 -0.92 \n cor.inc_lhd.rhd 8.860 2.977 0.95 -0.97 0.99 #not optimal"## check rePCA output =================================
summary(rePCA(m0_FL2))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5] [,6] [,7]
## Standard deviation 1.4261 0.8921 0.8329 0.52991 0.28957 0.001824 2.546e-05
## Proportion of Variance 0.5231 0.2047 0.1784 0.07222 0.02157 0.000000 0.000e+00
## Cumulative Proportion 0.5231 0.7278 0.9062 0.97843 1.00000 1.000000 1.000e+00
##
## $item
## Importance of components:
## [,1] [,2] [,3] [,4]
## Standard deviation 0.9041 0.16070 0.0003666 6.829e-05
## Proportion of Variance 0.9694 0.03063 0.0000000 0.000e+00
## Cumulative Proportion 0.9694 1.00000 1.0000000 1.000e+00# 5 of 7 components for Subj | 2 of 4 components for item contribute (variance)
# try model reduction of randEff ---------------------------------
## start with ZERO CORRELATION model =================================
summary(zcp_FL2 <- lmer(s8reln2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt||Subj)+
(1+lhd.rhd+cor.inc+
cor.inc_lhd.rhd||item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## s8reln2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + cor.inc_rd.rpt + grpd.ungr_rd.rpt ||
## Subj) + (1 + lhd.rhd + cor.inc + cor.inc_lhd.rhd || item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8373.4 8502.4 -4161.7 8323.4 1264
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2302 -0.5851 0.0076 0.5721 4.0602
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subj (Intercept) 1.596e+01 3.994e+00
## Subj.1 cor.inc 0.000e+00 0.000e+00
## Subj.2 grpd.ungr 7.066e+00 2.658e+00
## Subj.3 rd.rpt 1.214e+01 3.485e+00
## Subj.4 grpd.ungr_cor.inc 6.061e+01 7.785e+00
## Subj.5 cor.inc_rd.rpt 9.206e+00 3.034e+00
## Subj.6 grpd.ungr_rd.rpt 1.050e+01 3.240e+00
## item (Intercept) 9.079e+00 3.013e+00
## item.1 lhd.rhd 2.419e+00 1.555e+00
## item.2 cor.inc 1.152e+00 1.073e+00
## item.3 cor.inc_lhd.rhd 3.534e-09 5.945e-05
## Residual 2.858e+01 5.346e+00
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 43.31334 0.97786 53.26175 44.294 < 2e-16 ***
## grpd.ungr 4.25572 1.39282 31.88133 3.055 0.00452 **
## lhd.rhd 3.91615 1.55288 35.95425 2.522 0.01624 *
## cor.inc -0.01469 0.54154 54.67901 -0.027 0.97845
## rd.rpt 1.43481 1.48291 37.65768 0.968 0.33944
## grpd.ungr_lhd.rhd -1.10350 1.42861 35.65921 -0.772 0.44495
## grpd.ungr_cor.inc 10.09133 1.72258 33.76543 5.858 1.35e-06 ***
## cor.inc_lhd.rhd -1.56023 0.98699 443.50904 -1.581 0.11464
## grpd.ungr_rd.rpt -0.81302 2.65111 26.98838 -0.307 0.76145
## lhd.rhd_rd.rpt -4.54409 1.77172 36.27850 -2.565 0.01460 *
## cor.inc_rd.rpt 0.01200 1.17695 27.03668 0.010 0.99194
## grpd.ungr_lhd.rhd_rd.rpt -0.67822 2.31076 28.99594 -0.294 0.77123
## cor.inc_lhd.rhd_rd.rpt 4.29135 2.18515 31.85738 1.964 0.05832 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
'Random effects:
Groups Name Variance Std.Dev.
Subj (Intercept) 1.596e+01 3.9944841
Subj.1 cor.inc 0.000e+00 0.0000000 #no variance
Subj.2 grpd.ungr 7.065e+00 2.6579804
Subj.3 rd.rpt 1.214e+01 3.4847284
Subj.4 grpd.ungr_cor.inc 6.059e+01 7.7837605
Subj.5 cor.inc_rd.rpt 9.206e+00 3.0340750
Subj.6 grpd.ungr_rd.rpt 1.050e+01 3.2398369
item (Intercept) 9.081e+00 3.0134095
item.1 lhd.rhd 2.418e+00 1.5549853
item.2 cor.inc 1.152e+00 1.0732468
item.3 cor.inc_lhd.rhd 3.940e-08 0.0001985 #no variance
Residual 2.858e+01 5.3458264
Number of obs: 1289, groups: Subj, 32; item, 24'## [1] "Random effects:\n Groups Name Variance Std.Dev. \n Subj (Intercept) 1.596e+01 3.9944841\n Subj.1 cor.inc 0.000e+00 0.0000000 #no variance\n Subj.2 grpd.ungr 7.065e+00 2.6579804\n Subj.3 rd.rpt 1.214e+01 3.4847284\n Subj.4 grpd.ungr_cor.inc 6.059e+01 7.7837605\n Subj.5 cor.inc_rd.rpt 9.206e+00 3.0340750\n Subj.6 grpd.ungr_rd.rpt 1.050e+01 3.2398369\n item (Intercept) 9.081e+00 3.0134095\n item.1 lhd.rhd 2.418e+00 1.5549853\n item.2 cor.inc 1.152e+00 1.0732468\n item.3 cor.inc_lhd.rhd 3.940e-08 0.0001985 #no variance\n Residual 2.858e+01 5.3458264\n Number of obs: 1289, groups: Subj, 32; item, 24"## ANOVA model comparison - initial and zcp model =================================
anova(m0_FL2, zcp_FL2)# m0_FL2 has signif. better fit - removing all correlations too harsh
## back to initial model =================================
# try exclusion of randEff slopes with regards to zcp + initial model output
### FIRST #############################
# try exclusion of "cor.inc" as slope per Subj
summary(m1_FL2 <- lmer(s8reln2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd+cor.inc+
cor.inc_lhd.rhd|item),
REML = FALSE, pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## s8reln2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + grpd.ungr + rd.rpt + grpd.ungr_cor.inc +
## cor.inc_rd.rpt + grpd.ungr_rd.rpt | Subj) + (1 + lhd.rhd +
## cor.inc + cor.inc_lhd.rhd | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8352.4 8584.6 -4131.2 8262.4 1244
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2292 -0.6067 -0.0094 0.5770 4.0329
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 16.212 4.026
## grpd.ungr 7.284 2.699 -0.73
## rd.rpt 11.654 3.414 -0.01 0.05
## grpd.ungr_cor.inc 59.585 7.719 -0.19 0.08 -0.54
## cor.inc_rd.rpt 10.642 3.262 0.12 -0.01 -0.09 0.20
## grpd.ungr_rd.rpt 10.048 3.170 -0.28 -0.10 -0.11 0.21 -0.87
## item (Intercept) 6.724 2.593
## lhd.rhd 4.784 2.187 -1.00
## cor.inc 2.498 1.580 0.93 -0.93
## cor.inc_lhd.rhd 8.017 2.831 0.98 -0.98 0.99
## Residual 28.009 5.292
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 43.23713 0.92665 51.23790 46.659 < 2e-16 ***
## grpd.ungr 4.10716 1.20773 36.79100 3.401 0.00163 **
## lhd.rhd 3.84004 1.57169 37.67393 2.443 0.01936 *
## cor.inc 0.15422 0.55987 36.89006 0.275 0.78450
## rd.rpt 1.34602 1.32563 37.57666 1.015 0.31642
## grpd.ungr_lhd.rhd -0.75095 1.50785 41.50291 -0.498 0.62109
## grpd.ungr_cor.inc 10.15091 1.69708 26.76707 5.981 2.3e-06 ***
## cor.inc_lhd.rhd -1.66840 1.08237 52.68549 -1.541 0.12920
## grpd.ungr_rd.rpt 0.05174 2.19425 36.05503 0.024 0.98132
## lhd.rhd_rd.rpt -4.50130 1.71621 37.79469 -2.623 0.01250 *
## cor.inc_rd.rpt 0.38232 1.26105 29.98829 0.303 0.76385
## grpd.ungr_lhd.rhd_rd.rpt -1.50242 2.09139 35.17177 -0.718 0.47726
## cor.inc_lhd.rhd_rd.rpt 5.25565 2.44282 37.11923 2.151 0.03801 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
# Model failed to converge with 1 negative eigenvalue: -7.0e-01
### SECOND #############################
# back to initial model
# try exclusion of "cor.inc_lhd.rhd" as slope per item
summary(m2_FL2 <- lmer(s8reln2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd+cor.inc|item),
REML = FALSE, pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## s8reln2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + cor.inc_rd.rpt + grpd.ungr_rd.rpt | Subj) +
## (1 + lhd.rhd + cor.inc | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8365.0 8612.7 -4134.5 8269.0 1241
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2583 -0.6091 0.0013 0.5985 3.9704
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 17.448 4.177
## cor.inc 2.454 1.567 -0.43
## grpd.ungr 6.652 2.579 -0.75 0.21
## rd.rpt 10.095 3.177 0.02 0.47 -0.06
## grpd.ungr_cor.inc 49.346 7.025 -0.12 -0.84 0.23 -0.43
## cor.inc_rd.rpt 10.664 3.266 0.09 -0.07 0.18 -0.09 0.02
## grpd.ungr_rd.rpt 10.598 3.256 -0.16 -0.25 -0.20 -0.23 0.35
## item (Intercept) 7.427 2.725
## lhd.rhd 2.137 1.462 -0.98
## cor.inc 1.300 1.140 0.97 -0.92
## Residual 28.252 5.315
##
##
##
##
##
##
##
## -0.88
##
##
##
##
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 43.4684 0.9614 49.1543 45.215 < 2e-16 ***
## grpd.ungr 4.5384 1.2737 31.0033 3.563 0.00121 **
## lhd.rhd 3.8720 1.5907 32.2861 2.434 0.02063 *
## cor.inc -0.1135 0.5848 32.5070 -0.194 0.84736
## rd.rpt 1.6312 1.3516 35.7998 1.207 0.23539
## grpd.ungr_lhd.rhd -0.8135 1.3748 35.5770 -0.592 0.55775
## grpd.ungr_cor.inc 9.9653 1.5742 26.8695 6.330 9.1e-07 ***
## cor.inc_lhd.rhd -1.5300 1.0086 54.7002 -1.517 0.13502
## grpd.ungr_rd.rpt -0.2401 2.3384 30.5457 -0.103 0.91889
## lhd.rhd_rd.rpt -4.3291 1.5665 32.0087 -2.764 0.00940 **
## cor.inc_rd.rpt -0.4748 1.1800 30.6689 -0.402 0.69018
## grpd.ungr_lhd.rhd_rd.rpt -1.8997 2.1195 29.9436 -0.896 0.37725
## cor.inc_lhd.rhd_rd.rpt 4.7796 2.1681 31.2683 2.204 0.03498 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
# Model failed to converge with 1 negative eigenvalue: -4.8e-03
### THIRD #############################
# back to initial model
# try exclude "lhd.rhd" as slope for item - had corr of -1 in "m0_FL2"
summary(m3_FL2 <- lmer(s8reln2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt|Subj)+
(1+cor.inc+
cor.inc_lhd.rhd|item),
REML = FALSE, pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## s8reln2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + cor.inc_rd.rpt + grpd.ungr_rd.rpt | Subj) +
## (1 + cor.inc + cor.inc_lhd.rhd | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8384.1 8631.9 -4144.1 8288.1 1241
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3323 -0.5920 0.0206 0.5999 4.0344
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 17.3024 4.1596
## cor.inc 2.5649 1.6015 -0.38
## grpd.ungr 6.7412 2.5964 -0.76 0.14
## rd.rpt 9.9679 3.1572 0.00 0.43 -0.05
## grpd.ungr_cor.inc 45.5727 6.7508 -0.12 -0.86 0.24 -0.39
## cor.inc_rd.rpt 10.4548 3.2334 0.15 -0.08 0.14 -0.07 -0.01
## grpd.ungr_rd.rpt 10.6163 3.2583 -0.13 -0.28 -0.20 -0.25 0.38
## item (Intercept) 8.4238 2.9024
## cor.inc 1.5891 1.2606 0.76
## cor.inc_lhd.rhd 0.9492 0.9743 0.39 0.90
## Residual 28.6977 5.3570
##
##
##
##
##
##
##
## -0.88
##
##
##
##
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 43.4922 0.9811 49.8430 44.330 < 2e-16 ***
## grpd.ungr 4.4983 1.3402 30.4766 3.357 0.00213 **
## lhd.rhd 3.9094 1.5588 30.2495 2.508 0.01773 *
## cor.inc -0.1059 0.6006 28.2418 -0.176 0.86133
## rd.rpt 1.5621 1.4117 34.7016 1.107 0.27611
## grpd.ungr_lhd.rhd -0.8444 1.2487 29.7290 -0.676 0.50412
## grpd.ungr_cor.inc 10.2051 1.5421 26.0626 6.618 5.04e-07 ***
## cor.inc_lhd.rhd -1.3849 1.0358 44.8828 -1.337 0.18796
## grpd.ungr_rd.rpt -0.4127 2.5300 28.2477 -0.163 0.87158
## lhd.rhd_rd.rpt -4.4772 1.4662 26.0444 -3.054 0.00516 **
## cor.inc_rd.rpt -0.5320 1.2029 27.4099 -0.442 0.66177
## grpd.ungr_lhd.rhd_rd.rpt -1.3452 1.8993 32.8101 -0.708 0.48379
## cor.inc_lhd.rhd_rd.rpt 4.7692 2.2107 29.2855 2.157 0.03932 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular / Model failed to converge with 1 negative eigenvalue: -1.3e-03
# FINAL model ---------------------------------
## use initial model as final model =================================
# all model reduction steps lead to convergence fails
## initial model with factors and * for fixed effects =================================
# use REML=TRUE
summary(FL2 <- lmer(s8reln2 ~
condition*accuracy+condition*group*exp+
group*exp*accuracy+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd+cor.inc+
cor.inc_lhd.rhd|item),
REML = TRUE, pDat_mm))## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: s8reln2 ~ condition * accuracy + condition * group * exp + group *
## exp * accuracy + (1 + cor.inc + grpd.ungr + rd.rpt + grpd.ungr_cor.inc +
## cor.inc_rd.rpt + grpd.ungr_rd.rpt | Subj) + (1 + lhd.rhd +
## cor.inc + cor.inc_lhd.rhd | item)
## Data: pDat_mm
##
## REML criterion at convergence: 8222.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2224 -0.6160 -0.0031 0.5752 4.0520
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 18.900 4.347
## cor.inc 2.964 1.722 -0.47
## grpd.ungr 7.277 2.698 -0.72 0.29
## rd.rpt 10.781 3.283 0.05 0.47 -0.07
## grpd.ungr_cor.inc 51.359 7.167 -0.13 -0.80 0.18 -0.43
## cor.inc_rd.rpt 13.132 3.624 -0.02 -0.02 0.18 -0.10 0.05
## grpd.ungr_rd.rpt 12.222 3.496 -0.13 -0.31 -0.19 -0.24 0.39
## item (Intercept) 7.814 2.795
## lhd.rhd 5.639 2.375 -0.99
## cor.inc 3.376 1.837 0.83 -0.88
## cor.inc_lhd.rhd 10.222 3.197 0.92 -0.96 0.98
## Residual 27.752 5.268
##
##
##
##
##
##
##
## -0.84
##
##
##
##
##
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 43.4667 0.9946 45.0080 43.701 < 2e-16 ***
## condition1 4.7002 1.2599 31.2660 3.731 0.00076 ***
## accuracy1 -0.1078 0.6757 28.6669 -0.160 0.87434
## group1 3.9782 1.6939 32.7308 2.349 0.02504 *
## exp1 1.7214 1.3872 31.4757 1.241 0.22378
## condition1:accuracy1 9.5091 1.6415 25.8839 5.793 4.28e-06 ***
## condition1:group1 -0.4451 1.5410 38.1461 -0.289 0.77426
## condition1:exp1 0.1995 2.2935 28.0783 0.087 0.93132
## group1:exp1 -4.2525 1.7861 34.6633 -2.381 0.02290 *
## accuracy1:group1 -1.7256 1.2535 34.8840 -1.377 0.17738
## accuracy1:exp1 -0.7481 1.3588 22.3199 -0.551 0.58740
## condition1:group1:exp1 -2.2163 2.1720 30.4381 -1.020 0.31558
## accuracy1:group1:exp1 4.3406 2.6122 24.0041 1.662 0.10959
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
save(FL2, file = "FL2.Rdata")
# *Nice table output* ---------------------------------
# sjPlot::tab_model((FL2), show.se = TRUE, show.stat = TRUE, collapse.ci=TRUE, file="FL2.doc")
# clean workspace ---------------------------------
# ls()
rm(m0_FL2, m1_FL2, m2_FL2, m3_FL2, zcp_FL2)4.3 Model reduction - f0-range (F02)
# all "REML = FALSE" to use ANOVA model comparison later on
# INITIAL model with indicator variables ---------------------------------
summary(m0_F02 <- lmer(f0_range2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd+cor.inc+
cor.inc_lhd.rhd|item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## f0_range2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + cor.inc_rd.rpt + grpd.ungr_rd.rpt | Subj) +
## (1 + lhd.rhd + cor.inc + cor.inc_lhd.rhd | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 5955.8 6222.2 -2925.9 5851.8 1189
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4364 -0.4503 -0.0323 0.3961 6.8894
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 3.05305 1.7473
## cor.inc 0.06531 0.2556 1.00
## grpd.ungr 1.35506 1.1641 0.48 0.48
## rd.rpt 1.75924 1.3264 -0.12 -0.12 0.09
## grpd.ungr_cor.inc 7.10784 2.6661 0.49 0.49 -0.13 -0.37
## cor.inc_rd.rpt 0.93138 0.9651 0.83 0.83 0.55 -0.20 0.15
## grpd.ungr_rd.rpt 2.31544 1.5217 0.01 0.01 0.60 0.81 -0.57
## item (Intercept) 0.03354 0.1831
## lhd.rhd 0.52521 0.7247 -0.86
## cor.inc 0.12813 0.3580 0.39 -0.81
## cor.inc_lhd.rhd 0.40600 0.6372 0.92 -0.59 0.01
## Residual 5.43790 2.3319
##
##
##
##
##
##
##
## 0.09
##
##
##
##
##
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.3222 0.3313 32.2559 16.064 < 2e-16 ***
## grpd.ungr 2.0991 0.2890 28.6121 7.263 5.79e-08 ***
## lhd.rhd 1.1229 0.6330 37.7712 1.774 0.084135 .
## cor.inc 0.3210 0.2230 77.4360 1.440 0.154028
## rd.rpt 2.3654 0.3295 29.4105 7.178 6.16e-08 ***
## grpd.ungr_lhd.rhd 0.2354 0.6014 35.0779 0.391 0.697878
## grpd.ungr_cor.inc 3.0438 0.6349 26.6829 4.794 5.45e-05 ***
## cor.inc_lhd.rhd -0.3729 0.4349 55.4545 -0.858 0.394835
## grpd.ungr_rd.rpt 0.9006 0.4624 37.1558 1.948 0.059030 .
## lhd.rhd_rd.rpt 2.9512 0.6819 36.1888 4.328 0.000114 ***
## cor.inc_rd.rpt 0.2083 0.4645 27.9049 0.448 0.657271
## grpd.ungr_lhd.rhd_rd.rpt 1.8301 0.9141 39.4978 2.002 0.052184 .
## cor.inc_lhd.rhd_rd.rpt 1.2205 0.9171 23.6050 1.331 0.195959
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
# Model failed to converge with 2 negative eigenvalues: -1.5e-03 -1.6e+01
summary(rePCA(m0_F02))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5] [,6] [,7]
## Standard deviation 1.3133 0.8990 0.6007 0.3659 0.15069 4.165e-05 1.417e-17
## Proportion of Variance 0.5654 0.2650 0.1183 0.0439 0.00744 0.000e+00 0.000e+00
## Cumulative Proportion 0.5654 0.8304 0.9487 0.9926 1.00000 1.000e+00 1.000e+00
##
## $item
## Importance of components:
## [,1] [,2] [,3] [,4]
## Standard deviation 0.3884 0.2238 0.0002056 0
## Proportion of Variance 0.7508 0.2492 0.0000000 0
## Cumulative Proportion 0.7508 1.0000 1.0000000 1# 5 of 7 components for Subj | 2 of 4 components for item contribute (variance)
# try model reduction of randEff ---------------------------------
## start with ZERO CORRELATION model =================================
summary(zcp_F02 <- lmer(f0_range2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt||Subj)+
(1+lhd.rhd+cor.inc+
cor.inc_lhd.rhd||item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## f0_range2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + cor.inc_rd.rpt + grpd.ungr_rd.rpt ||
## Subj) + (1 + lhd.rhd + cor.inc + cor.inc_lhd.rhd || item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 5949.1 6077.2 -2949.6 5899.1 1216
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.0695 -0.4312 -0.0440 0.3984 6.7809
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subj (Intercept) 3.238e+00 1.799e+00
## Subj.1 cor.inc 0.000e+00 0.000e+00
## Subj.2 grpd.ungr 1.220e+00 1.104e+00
## Subj.3 rd.rpt 1.469e+00 1.212e+00
## Subj.4 grpd.ungr_cor.inc 5.502e+00 2.346e+00
## Subj.5 cor.inc_rd.rpt 1.047e+00 1.023e+00
## Subj.6 grpd.ungr_rd.rpt 4.089e-01 6.394e-01
## item (Intercept) 4.275e-02 2.067e-01
## item.1 lhd.rhd 4.676e-01 6.838e-01
## item.2 cor.inc 0.000e+00 0.000e+00
## item.3 cor.inc_lhd.rhd 6.944e-09 8.333e-05
## Residual 5.631e+00 2.373e+00
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.3717 0.3417 35.2945 15.722 < 2e-16 ***
## grpd.ungr 2.0899 0.2842 31.3584 7.354 2.61e-08 ***
## lhd.rhd 1.6858 0.6924 36.9313 2.435 0.019847 *
## cor.inc 0.4174 0.2115 411.0690 1.974 0.049099 *
## rd.rpt 2.2444 0.3228 36.6579 6.954 3.43e-08 ***
## grpd.ungr_lhd.rhd 0.1203 0.5973 33.5177 0.201 0.841639
## grpd.ungr_cor.inc 3.2277 0.5741 32.5001 5.622 3.09e-06 ***
## cor.inc_lhd.rhd -0.2549 0.4215 430.0617 -0.605 0.545551
## grpd.ungr_rd.rpt 0.7294 0.3836 23.4116 1.902 0.069589 .
## lhd.rhd_rd.rpt 2.5498 0.6815 39.9618 3.741 0.000575 ***
## cor.inc_rd.rpt 0.2252 0.4422 29.2211 0.509 0.614403
## grpd.ungr_lhd.rhd_rd.rpt 1.3568 0.8804 23.7869 1.541 0.136478
## cor.inc_lhd.rhd_rd.rpt 1.2258 0.8847 29.2133 1.386 0.176386
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
summary(zcp_F02)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## f0_range2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + cor.inc_rd.rpt + grpd.ungr_rd.rpt ||
## Subj) + (1 + lhd.rhd + cor.inc + cor.inc_lhd.rhd || item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 5949.1 6077.2 -2949.6 5899.1 1216
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.0695 -0.4312 -0.0440 0.3984 6.7809
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subj (Intercept) 3.238e+00 1.799e+00
## Subj.1 cor.inc 0.000e+00 0.000e+00
## Subj.2 grpd.ungr 1.220e+00 1.104e+00
## Subj.3 rd.rpt 1.469e+00 1.212e+00
## Subj.4 grpd.ungr_cor.inc 5.502e+00 2.346e+00
## Subj.5 cor.inc_rd.rpt 1.047e+00 1.023e+00
## Subj.6 grpd.ungr_rd.rpt 4.089e-01 6.394e-01
## item (Intercept) 4.275e-02 2.067e-01
## item.1 lhd.rhd 4.676e-01 6.838e-01
## item.2 cor.inc 0.000e+00 0.000e+00
## item.3 cor.inc_lhd.rhd 6.944e-09 8.333e-05
## Residual 5.631e+00 2.373e+00
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.3717 0.3417 35.2945 15.722 < 2e-16 ***
## grpd.ungr 2.0899 0.2842 31.3584 7.354 2.61e-08 ***
## lhd.rhd 1.6858 0.6924 36.9313 2.435 0.019847 *
## cor.inc 0.4174 0.2115 411.0690 1.974 0.049099 *
## rd.rpt 2.2444 0.3228 36.6579 6.954 3.43e-08 ***
## grpd.ungr_lhd.rhd 0.1203 0.5973 33.5177 0.201 0.841639
## grpd.ungr_cor.inc 3.2277 0.5741 32.5001 5.622 3.09e-06 ***
## cor.inc_lhd.rhd -0.2549 0.4215 430.0617 -0.605 0.545551
## grpd.ungr_rd.rpt 0.7294 0.3836 23.4116 1.902 0.069589 .
## lhd.rhd_rd.rpt 2.5498 0.6815 39.9618 3.741 0.000575 ***
## cor.inc_rd.rpt 0.2252 0.4422 29.2211 0.509 0.614403
## grpd.ungr_lhd.rhd_rd.rpt 1.3568 0.8804 23.7869 1.541 0.136478
## cor.inc_lhd.rhd_rd.rpt 1.2258 0.8847 29.2133 1.386 0.176386
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular') 'Random effects:
Groups Name Variance Std.Dev.
Subj (Intercept) 3.238e+00 1.7993344
Subj.1 cor.inc 0.000e+00 0.0000000 # no variance
Subj.2 grpd.ungr 1.220e+00 1.1043709
Subj.3 rd.rpt 1.469e+00 1.2118683
Subj.4 grpd.ungr_cor.inc 5.502e+00 2.3457306
Subj.5 cor.inc_rd.rpt 1.047e+00 1.0231915 # mean var. + sd almost equal
Subj.6 grpd.ungr_rd.rpt 4.087e-01 0.6392808 # mean var. smaller than sd
item (Intercept) 4.275e-02 0.2067492
item.1 lhd.rhd 4.676e-01 0.6838274
item.2 cor.inc 8.172e-08 0.0002859 # no variance
item.3 cor.inc_lhd.rhd 0.000e+00 0.0000000 # no variance
Residual 5.631e+00 2.3728895
Number of obs: 1241, groups: Subj, 32; item, 24'## [1] "Random effects:\n Groups Name Variance Std.Dev. \n Subj (Intercept) 3.238e+00 1.7993344\n Subj.1 cor.inc 0.000e+00 0.0000000 # no variance\n Subj.2 grpd.ungr 1.220e+00 1.1043709\n Subj.3 rd.rpt 1.469e+00 1.2118683\n Subj.4 grpd.ungr_cor.inc 5.502e+00 2.3457306\n Subj.5 cor.inc_rd.rpt 1.047e+00 1.0231915 # mean var. + sd almost equal\n Subj.6 grpd.ungr_rd.rpt 4.087e-01 0.6392808 # mean var. smaller than sd\n item (Intercept) 4.275e-02 0.2067492\n item.1 lhd.rhd 4.676e-01 0.6838274\n item.2 cor.inc 8.172e-08 0.0002859 # no variance\n item.3 cor.inc_lhd.rhd 0.000e+00 0.0000000 # no variance\n Residual 5.631e+00 2.3728895\n Number of obs: 1241, groups: Subj, 32; item, 24"## ANOVA model comparison - initial - zcp model =================================
anova(m0_F02, zcp_F02)# initial model m0_F02 significantly better variance resolution
# BUT zcp_F02 has smaller AIC - therefore better model fit
# reasonable to look further into zcp model
## iterative reduction of randEff =================================
### FIRST #############################
# try exclusion of "cor.inc" as slope per Subj
summary(zcp1_F02 <- lmer(f0_range2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt||Subj)+
(1+lhd.rhd+cor.inc+
cor.inc_lhd.rhd||item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## f0_range2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + grpd.ungr + rd.rpt + grpd.ungr_cor.inc +
## cor.inc_rd.rpt + grpd.ungr_rd.rpt || Subj) + (1 + lhd.rhd +
## cor.inc + cor.inc_lhd.rhd || item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 5947.1 6070.1 -2949.6 5899.1 1217
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.0695 -0.4312 -0.0440 0.3985 6.7809
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subj (Intercept) 3.23785 1.7994
## Subj.1 grpd.ungr 1.21953 1.1043
## Subj.2 rd.rpt 1.46871 1.2119
## Subj.3 grpd.ungr_cor.inc 5.50201 2.3456
## Subj.4 cor.inc_rd.rpt 1.04687 1.0232
## Subj.5 grpd.ungr_rd.rpt 0.40876 0.6393
## item (Intercept) 0.04275 0.2068
## item.1 lhd.rhd 0.46757 0.6838
## item.2 cor.inc 0.00000 0.0000
## item.3 cor.inc_lhd.rhd 0.00000 0.0000
## Residual 5.63060 2.3729
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.3717 0.3417 35.2937 15.722 < 2e-16 ***
## grpd.ungr 2.0899 0.2842 31.3596 7.354 2.61e-08 ***
## lhd.rhd 1.6858 0.6924 36.9304 2.435 0.019848 *
## cor.inc 0.4174 0.2115 411.0560 1.974 0.049099 *
## rd.rpt 2.2444 0.3228 36.6586 6.954 3.43e-08 ***
## grpd.ungr_lhd.rhd 0.1203 0.5973 33.5181 0.201 0.841642
## grpd.ungr_cor.inc 3.2277 0.5741 32.5008 5.622 3.09e-06 ***
## cor.inc_lhd.rhd -0.2550 0.4214 430.0484 -0.605 0.545537
## grpd.ungr_rd.rpt 0.7294 0.3836 23.4103 1.902 0.069589 .
## lhd.rhd_rd.rpt 2.5498 0.6815 39.9620 3.741 0.000575 ***
## cor.inc_rd.rpt 0.2252 0.4422 29.2201 0.509 0.614390
## grpd.ungr_lhd.rhd_rd.rpt 1.3569 0.8804 23.7853 1.541 0.136469
## cor.inc_lhd.rhd_rd.rpt 1.2258 0.8847 29.2124 1.386 0.176374
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
'Random effects:
Groups Name Variance Std.Dev.
Subj (Intercept) 3.23828 1.7995
Subj.1 grpd.ungr 1.21951 1.1043
Subj.2 rd.rpt 1.46867 1.2119
Subj.3 grpd.ungr_cor.inc 5.50127 2.3455
Subj.4 cor.inc_rd.rpt 1.04691 1.0232 # mean var. + sd almost equal
Subj.5 grpd.ungr_rd.rpt 0.40879 0.6394 # mean var. smaller than sd
item (Intercept) 0.04276 0.2068
item.1 lhd.rhd 0.46757 0.6838
item.2 cor.inc 0.00000 0.0000 # no variance
item.3 cor.inc_lhd.rhd 0.00000 0.0000 # no variance
Residual 5.63059 2.3729
Number of obs: 1241, groups: Subj, 32; item, 24'## [1] "Random effects:\n Groups Name Variance Std.Dev.\n Subj (Intercept) 3.23828 1.7995 \n Subj.1 grpd.ungr 1.21951 1.1043 \n Subj.2 rd.rpt 1.46867 1.2119 \n Subj.3 grpd.ungr_cor.inc 5.50127 2.3455 \n Subj.4 cor.inc_rd.rpt 1.04691 1.0232 # mean var. + sd almost equal \n Subj.5 grpd.ungr_rd.rpt 0.40879 0.6394 # mean var. smaller than sd \n item (Intercept) 0.04276 0.2068 \n item.1 lhd.rhd 0.46757 0.6838 \n item.2 cor.inc 0.00000 0.0000 # no variance \n item.3 cor.inc_lhd.rhd 0.00000 0.0000 # no variance\n Residual 5.63059 2.3729 \n Number of obs: 1241, groups: Subj, 32; item, 24"summary(rePCA(zcp1_F02))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5] [,6]
## Standard deviation 0.9885 0.7583 0.5107 0.46539 0.43119 0.26944
## Proportion of Variance 0.4270 0.2513 0.1140 0.09466 0.08125 0.03173
## Cumulative Proportion 0.4270 0.6784 0.7924 0.88702 0.96827 1.00000
##
## $item
## Importance of components:
## [,1] [,2] [,3] [,4]
## Standard deviation 0.2882 0.08714 0 0
## Proportion of Variance 0.9162 0.08378 0 0
## Cumulative Proportion 0.9162 1.00000 1 1# components for Subj seem okay | 2 of 4 components for item have zero variance
### ANOVA model comparison - reduced zcp - zcp model #############################
anova(zcp_F02, zcp1_F02)# not significant, AIC smaller for zcp1_F02 - preferred
### SECOND #############################
# try also exclusion also of "cor.inc_lhd.rhd" as slope per item
summary(zcp2_F02 <- lmer(f0_range2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt||Subj)+
(1+lhd.rhd+cor.inc||item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## f0_range2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + grpd.ungr + rd.rpt + grpd.ungr_cor.inc +
## cor.inc_rd.rpt + grpd.ungr_rd.rpt || Subj) + (1 + lhd.rhd +
## cor.inc || item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 5945.1 6062.9 -2949.6 5899.1 1218
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.0695 -0.4312 -0.0440 0.3984 6.7809
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subj (Intercept) 3.238e+00 1.799e+00
## Subj.1 grpd.ungr 1.219e+00 1.104e+00
## Subj.2 rd.rpt 1.469e+00 1.212e+00
## Subj.3 grpd.ungr_cor.inc 5.503e+00 2.346e+00
## Subj.4 cor.inc_rd.rpt 1.047e+00 1.023e+00
## Subj.5 grpd.ungr_rd.rpt 4.090e-01 6.395e-01
## item (Intercept) 4.274e-02 2.067e-01
## item.1 lhd.rhd 4.675e-01 6.838e-01
## item.2 cor.inc 4.418e-09 6.647e-05
## Residual 5.631e+00 2.373e+00
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.3717 0.3417 35.2942 15.722 < 2e-16 ***
## grpd.ungr 2.0899 0.2842 31.3593 7.354 2.61e-08 ***
## lhd.rhd 1.6858 0.6924 36.9310 2.435 0.019847 *
## cor.inc 0.4174 0.2115 411.0646 1.974 0.049091 *
## rd.rpt 2.2444 0.3228 36.6570 6.953 3.43e-08 ***
## grpd.ungr_lhd.rhd 0.1203 0.5973 33.5188 0.201 0.841652
## grpd.ungr_cor.inc 3.2277 0.5741 32.4985 5.622 3.09e-06 ***
## cor.inc_lhd.rhd -0.2549 0.4215 430.0566 -0.605 0.545597
## grpd.ungr_rd.rpt 0.7294 0.3836 23.4124 1.902 0.069593 .
## lhd.rhd_rd.rpt 2.5498 0.6815 39.9611 3.741 0.000575 ***
## cor.inc_rd.rpt 0.2252 0.4422 29.2213 0.509 0.614389
## grpd.ungr_lhd.rhd_rd.rpt 1.3568 0.8804 23.7880 1.541 0.136493
## cor.inc_lhd.rhd_rd.rpt 1.2258 0.8847 29.2135 1.386 0.176380
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
'Random effects:
Groups Name Variance Std.Dev.
Subj (Intercept) 3.23777 1.7994
Subj.1 grpd.ungr 1.21950 1.1043
Subj.2 rd.rpt 1.46869 1.2119
Subj.3 grpd.ungr_cor.inc 5.50204 2.3456
Subj.4 cor.inc_rd.rpt 1.04707 1.0233 # mean var. + sd almost equal
Subj.5 grpd.ungr_rd.rpt 0.40879 0.6394 # mean var. smaller than sd
item (Intercept) 0.04275 0.2068
item.1 lhd.rhd 0.46755 0.6838
item.2 cor.inc 0.00000 0.0000 # no variance
Residual 5.63061 2.3729
Number of obs: 1241, groups: Subj, 32; item, 24'## [1] "Random effects:\n Groups Name Variance Std.Dev.\n Subj (Intercept) 3.23777 1.7994 \n Subj.1 grpd.ungr 1.21950 1.1043 \n Subj.2 rd.rpt 1.46869 1.2119 \n Subj.3 grpd.ungr_cor.inc 5.50204 2.3456 \n Subj.4 cor.inc_rd.rpt 1.04707 1.0233 # mean var. + sd almost equal \n Subj.5 grpd.ungr_rd.rpt 0.40879 0.6394 # mean var. smaller than sd\n item (Intercept) 0.04275 0.2068 \n item.1 lhd.rhd 0.46755 0.6838 \n item.2 cor.inc 0.00000 0.0000 # no variance\n Residual 5.63061 2.3729 \n Number of obs: 1241, groups: Subj, 32; item, 24"summary(rePCA(zcp2_F02))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5] [,6]
## Standard deviation 0.9886 0.7583 0.5107 0.46538 0.43127 0.26951
## Proportion of Variance 0.4271 0.2513 0.1140 0.09464 0.08128 0.03174
## Cumulative Proportion 0.4271 0.6784 0.7923 0.88698 0.96826 1.00000
##
## $item
## Importance of components:
## [,1] [,2] [,3]
## Standard deviation 0.2882 0.08713 2.801e-05
## Proportion of Variance 0.9162 0.08376 0.000e+00
## Cumulative Proportion 0.9162 1.00000 1.000e+00# components for Subj seem okay - 1 of 3 components for item has zero variance
### ANOVA model comparison - reduced 2nd zcp - reduced 1st zcp model #############################
anova(zcp2_F02, zcp1_F02)# not significant and AIC smaller for zcp2_F02 - preferred
### THIRD #############################
# try exclusion also of "cor.inc" as slope per item
summary(zcp3_F02 <- lmer(f0_range2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt||Subj)+
(1+lhd.rhd||item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## f0_range2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + grpd.ungr + rd.rpt + grpd.ungr_cor.inc +
## cor.inc_rd.rpt + grpd.ungr_rd.rpt || Subj) + (1 + lhd.rhd || item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 5943.1 6055.8 -2949.6 5899.1 1219
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.0695 -0.4312 -0.0440 0.3984 6.7809
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subj (Intercept) 3.23798 1.7994
## Subj.1 grpd.ungr 1.21956 1.1043
## Subj.2 rd.rpt 1.46878 1.2119
## Subj.3 grpd.ungr_cor.inc 5.50252 2.3457
## Subj.4 cor.inc_rd.rpt 1.04718 1.0233
## Subj.5 grpd.ungr_rd.rpt 0.40874 0.6393
## item (Intercept) 0.04274 0.2067
## item.1 lhd.rhd 0.46754 0.6838
## Residual 5.63059 2.3729
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.3717 0.3417 35.2911 15.721 < 2e-16 ***
## grpd.ungr 2.0899 0.2842 31.3580 7.354 2.61e-08 ***
## lhd.rhd 1.6858 0.6924 36.9278 2.435 0.019851 *
## cor.inc 0.4174 0.2115 411.0610 1.974 0.049095 *
## rd.rpt 2.2444 0.3228 36.6557 6.953 3.43e-08 ***
## grpd.ungr_lhd.rhd 0.1203 0.5973 33.5178 0.201 0.841660
## grpd.ungr_cor.inc 3.2277 0.5741 32.4987 5.622 3.09e-06 ***
## cor.inc_lhd.rhd -0.2549 0.4215 430.0530 -0.605 0.545577
## grpd.ungr_rd.rpt 0.7294 0.3835 23.4097 1.902 0.069579 .
## lhd.rhd_rd.rpt 2.5498 0.6815 39.9601 3.741 0.000575 ***
## cor.inc_rd.rpt 0.2252 0.4422 29.2210 0.509 0.614373
## grpd.ungr_lhd.rhd_rd.rpt 1.3568 0.8804 23.7859 1.541 0.136476
## cor.inc_lhd.rhd_rd.rpt 1.2258 0.8847 29.2132 1.386 0.176381
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 'Random effects:
Groups Name Variance Std.Dev.
Subj (Intercept) 3.23788 1.7994
Subj.1 grpd.ungr 1.21961 1.1044
Subj.2 rd.rpt 1.46866 1.2119
Subj.3 grpd.ungr_cor.inc 5.50221 2.3457
Subj.4 cor.inc_rd.rpt 1.04710 1.0233 # mean var. + sd almost equal
Subj.5 grpd.ungr_rd.rpt 0.40869 0.6393 # mean var. smaller than sd
item (Intercept) 0.04274 0.2067
item.1 lhd.rhd 0.46756 0.6838
Residual 5.63060 2.3729
Number of obs: 1241, groups: Subj, 32; item, 24'## [1] "Random effects:\n Groups Name Variance Std.Dev.\n Subj (Intercept) 3.23788 1.7994 \n Subj.1 grpd.ungr 1.21961 1.1044 \n Subj.2 rd.rpt 1.46866 1.2119 \n Subj.3 grpd.ungr_cor.inc 5.50221 2.3457 \n Subj.4 cor.inc_rd.rpt 1.04710 1.0233 # mean var. + sd almost equal \n Subj.5 grpd.ungr_rd.rpt 0.40869 0.6393 # mean var. smaller than sd \n item (Intercept) 0.04274 0.2067 \n item.1 lhd.rhd 0.46756 0.6838 \n Residual 5.63060 2.3729 \n Number of obs: 1241, groups: Subj, 32; item, 24"summary(rePCA(zcp3_F02))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5] [,6]
## Standard deviation 0.9886 0.7583 0.5107 0.46540 0.43125 0.26943
## Proportion of Variance 0.4271 0.2513 0.1140 0.09465 0.08127 0.03172
## Cumulative Proportion 0.4271 0.6784 0.7923 0.88701 0.96828 1.00000
##
## $item
## Importance of components:
## [,1] [,2]
## Standard deviation 0.2882 0.08712
## Proportion of Variance 0.9162 0.08376
## Cumulative Proportion 0.9162 1.00000# all components seem okay
### ANOVA model comparison - reduced 3rd zcp - reduced 2nd zcp model #############################
anova(zcp3_F02, zcp2_F02)# not significant and AIC smaller for zcp3_F02 - preferred
### FOURTH #############################
# try exclusion also of "grpd.ungr_rd.rpt" as slope per Subj
summary(zcp4_F02 <- lmer(f0_range2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt||Subj)+
(1+lhd.rhd||item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## f0_range2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + grpd.ungr + rd.rpt + grpd.ungr_cor.inc +
## cor.inc_rd.rpt || Subj) + (1 + lhd.rhd || item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 5941.5 6049.0 -2949.7 5899.5 1220
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.0130 -0.4368 -0.0452 0.4072 6.7353
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subj (Intercept) 3.2344 1.7984
## Subj.1 grpd.ungr 1.2075 1.0989
## Subj.2 rd.rpt 1.4482 1.2034
## Subj.3 grpd.ungr_cor.inc 5.4742 2.3397
## Subj.4 cor.inc_rd.rpt 1.1115 1.0543
## item (Intercept) 0.0418 0.2045
## item.1 lhd.rhd 0.4660 0.6827
## Residual 5.6543 2.3779
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.3698 0.3414 35.3045 15.727 < 2e-16 ***
## grpd.ungr 2.0785 0.2824 31.2622 7.359 2.63e-08 ***
## lhd.rhd 1.6888 0.6920 36.9575 2.441 0.019573 *
## cor.inc 0.4239 0.2110 416.2915 2.009 0.045174 *
## rd.rpt 2.2337 0.3214 36.6877 6.949 3.46e-08 ***
## grpd.ungr_lhd.rhd 0.1079 0.5942 33.4156 0.182 0.857023
## grpd.ungr_cor.inc 3.2464 0.5716 32.5122 5.680 2.61e-06 ***
## cor.inc_lhd.rhd -0.2559 0.4206 433.2385 -0.608 0.543277
## grpd.ungr_rd.rpt 0.7355 0.3617 37.1966 2.034 0.049169 *
## lhd.rhd_rd.rpt 2.5438 0.6793 40.0039 3.745 0.000569 ***
## cor.inc_rd.rpt 0.2514 0.4421 28.8661 0.569 0.573949
## grpd.ungr_lhd.rhd_rd.rpt 1.3749 0.8438 27.7374 1.629 0.114535
## cor.inc_lhd.rhd_rd.rpt 1.2140 0.8844 28.8625 1.373 0.180424
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.00466327 (tol = 0.002, component 1) 'Random effects:
Groups Name Variance Std.Dev.
Subj (Intercept) 3.23610 1.7989
Subj.1 grpd.ungr 1.20739 1.0988
Subj.2 rd.rpt 1.44618 1.2026
Subj.3 grpd.ungr_cor.inc 5.47603 2.3401
Subj.4 cor.inc_rd.rpt 1.11054 1.0538 # mean var. + sd almost equal
item (Intercept) 0.04183 0.2045
item.1 lhd.rhd 0.46619 0.6828
Residual 5.65429 2.3779
Number of obs: 1241, groups: Subj, 32; item, 24'## [1] "Random effects:\n Groups Name Variance Std.Dev.\n Subj (Intercept) 3.23610 1.7989 \n Subj.1 grpd.ungr 1.20739 1.0988 \n Subj.2 rd.rpt 1.44618 1.2026 \n Subj.3 grpd.ungr_cor.inc 5.47603 2.3401 \n Subj.4 cor.inc_rd.rpt 1.11054 1.0538 # mean var. + sd almost equal \n item (Intercept) 0.04183 0.2045 \n item.1 lhd.rhd 0.46619 0.6828 \n Residual 5.65429 2.3779 \n Number of obs: 1241, groups: Subj, 32; item, 24"summary(rePCA(zcp4_F02))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5]
## Standard deviation 0.9839 0.7563 0.5061 0.46212 0.44337
## Proportion of Variance 0.4388 0.2592 0.1161 0.09679 0.08909
## Cumulative Proportion 0.4388 0.6980 0.8141 0.91091 1.00000
##
## $item
## Importance of components:
## [,1] [,2]
## Standard deviation 0.2871 0.08598
## Proportion of Variance 0.9177 0.08231
## Cumulative Proportion 0.9177 1.00000# all components seem okay
### ANOVA model comparison - reduced 4th zcp - reduced 3rd zcp model #############################
anova(zcp4_F02, zcp3_F02)# not significant and AIC smaller for zcp4_F02 - preferred
### FIFTH #############################
# try exclusion also of "cor.inc_rd.rpt" as slope per Subj
summary(zcp5_F02 <- lmer(f0_range2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc||Subj)+
(1+lhd.rhd||item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## f0_range2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + grpd.ungr + rd.rpt + grpd.ungr_cor.inc ||
## Subj) + (1 + lhd.rhd || item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 5940.3 6042.8 -2950.2 5900.3 1221
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.0144 -0.4412 -0.0484 0.4113 6.7447
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subj (Intercept) 3.29998 1.8166
## Subj.1 grpd.ungr 1.16990 1.0816
## Subj.2 rd.rpt 1.59736 1.2639
## Subj.3 grpd.ungr_cor.inc 5.49785 2.3447
## item (Intercept) 0.04001 0.2000
## item.1 lhd.rhd 0.45381 0.6737
## Residual 5.69228 2.3858
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.36753 0.34430 35.47415 15.590 < 2e-16 ***
## grpd.ungr 2.06464 0.27855 31.26358 7.412 2.27e-08 ***
## lhd.rhd 1.68259 0.69746 37.06045 2.412 0.020910 *
## cor.inc 0.40069 0.20595 925.92923 1.946 0.052006 .
## rd.rpt 2.28622 0.32548 36.38732 7.024 2.88e-08 ***
## grpd.ungr_lhd.rhd 0.09631 0.58684 33.40822 0.164 0.870634
## grpd.ungr_cor.inc 3.24632 0.57027 32.65268 5.693 2.47e-06 ***
## cor.inc_lhd.rhd -0.28707 0.41072 963.49607 -0.699 0.484764
## grpd.ungr_rd.rpt 0.75711 0.35575 36.10486 2.128 0.040213 *
## lhd.rhd_rd.rpt 2.59721 0.68657 39.77712 3.783 0.000511 ***
## cor.inc_rd.rpt 0.21811 0.39229 1126.74603 0.556 0.578332
## grpd.ungr_lhd.rhd_rd.rpt 1.51546 0.83233 27.13287 1.821 0.079694 .
## cor.inc_lhd.rhd_rd.rpt 1.23852 0.78475 1121.37900 1.578 0.114794
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 'Random effects:
Groups Name Variance Std.Dev.
Subj (Intercept) 3.29999 1.8166
Subj.1 grpd.ungr 1.16990 1.0816
Subj.2 rd.rpt 1.59736 1.2639
Subj.3 grpd.ungr_cor.inc 5.49785 2.3447
item (Intercept) 0.04001 0.2000
item.1 lhd.rhd 0.45381 0.6737
Residual 5.69227 2.3858
Number of obs: 1241, groups: Subj, 32; item, 24'## [1] "Random effects:\n Groups Name Variance Std.Dev.\n Subj (Intercept) 3.29999 1.8166 \n Subj.1 grpd.ungr 1.16990 1.0816 \n Subj.2 rd.rpt 1.59736 1.2639 \n Subj.3 grpd.ungr_cor.inc 5.49785 2.3447 \n item (Intercept) 0.04001 0.2000 \n item.1 lhd.rhd 0.45381 0.6737 \n Residual 5.69227 2.3858 \n Number of obs: 1241, groups: Subj, 32; item, 24"summary(rePCA(zcp5_F02))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4]
## Standard deviation 0.9828 0.7614 0.5297 0.4533
## Proportion of Variance 0.4754 0.2853 0.1381 0.1012
## Cumulative Proportion 0.4754 0.7607 0.8988 1.0000
##
## $item
## Importance of components:
## [,1] [,2]
## Standard deviation 0.2824 0.08383
## Proportion of Variance 0.9190 0.08101
## Cumulative Proportion 0.9190 1.00000# all components seem okay
### ANOVA model comparison - reduced 5th zcp - reduced 4th zcp model #############################
anova(zcp5_F02, zcp4_F02)# not significant and AIC smaller for zcp5_F02 - preferred
# I go with zcp5_F02 as final zcp model and check correlations now
## REDUCED FINAL model with correlations =================================
summary(corr5_F02 <- lmer(f0_range2 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc|Subj)+
(1+lhd.rhd|item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## f0_range2 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + grpd.ungr + rd.rpt + grpd.ungr_cor.inc |
## Subj) + (1 + lhd.rhd | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 5934.6 6072.9 -2940.3 5880.6 1214
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.0620 -0.4399 -0.0492 0.4312 6.7337
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 3.42803 1.8515
## grpd.ungr 1.11676 1.0568 0.57
## rd.rpt 1.66202 1.2892 0.05 -0.09
## grpd.ungr_cor.inc 5.26044 2.2936 0.57 0.14 -0.04
## item (Intercept) 0.04727 0.2174
## lhd.rhd 0.37060 0.6088 -1.00
## Residual 5.68054 2.3834
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.26272 0.34923 33.40963 15.069 < 2e-16 ***
## grpd.ungr 2.05936 0.27595 29.61640 7.463 2.78e-08 ***
## lhd.rhd 0.92840 0.64789 38.75062 1.433 0.159884
## cor.inc 0.39009 0.20503 793.88170 1.903 0.057458 .
## rd.rpt 2.37051 0.32570 37.53587 7.278 1.11e-08 ***
## grpd.ungr_lhd.rhd 0.04822 0.57323 32.83986 0.084 0.933469
## grpd.ungr_cor.inc 3.23058 0.56081 29.55348 5.761 2.88e-06 ***
## cor.inc_lhd.rhd -0.25374 0.40878 828.65554 -0.621 0.534951
## grpd.ungr_rd.rpt 0.79556 0.36019 46.68496 2.209 0.032139 *
## lhd.rhd_rd.rpt 2.82898 0.67214 39.57005 4.209 0.000143 ***
## cor.inc_rd.rpt 0.29084 0.38932 987.57499 0.747 0.455216
## grpd.ungr_lhd.rhd_rd.rpt 1.41032 0.79387 31.45708 1.777 0.085324 .
## cor.inc_lhd.rhd_rd.rpt 1.24264 0.77915 989.06389 1.595 0.111060
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
### ANOVA model comparison - reduced final model with correlation - initial model #############################
anova(corr5_F02, m0_F02)# not significant and AIC smaller for corr5_F02 - preferred
# initial model does not accumulate more information
# FINAL model ---------------------------------
## model with factors and * for fixed effects =================================
# use REML=TRUE
summary(F02 <- lmer(f0_range2 ~
condition*accuracy+condition*group*exp+group*exp*accuracy+
(1+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc|Subj)+
(1+lhd.rhd|item),
REML = TRUE,pDat_mm))## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: f0_range2 ~ condition * accuracy + condition * group * exp +
## group * exp * accuracy + (1 + grpd.ungr + rd.rpt + grpd.ungr_cor.inc |
## Subj) + (1 + lhd.rhd | item)
## Data: pDat_mm
##
## REML criterion at convergence: 5878.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.0371 -0.4386 -0.0492 0.4286 6.7220
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 3.65139 1.9109
## grpd.ungr 1.27163 1.1277 0.55
## rd.rpt 1.83021 1.3529 0.05 -0.08
## grpd.ungr_cor.inc 5.58866 2.3640 0.56 0.11 -0.05
## item (Intercept) 0.05585 0.2363
## lhd.rhd 0.43847 0.6622 -1.00
## Residual 5.69990 2.3874
## Number of obs: 1241, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 5.26412 0.36027 31.76563 14.612 1.19e-15 ***
## condition1 2.06366 0.28833 28.28438 7.157 8.17e-08 ***
## accuracy1 0.38055 0.20676 789.12333 1.841 0.066063 .
## group1 0.93450 0.67099 36.46869 1.393 0.172136
## exp1 2.37714 0.33821 34.85543 7.029 3.59e-08 ***
## condition1:accuracy1 3.20448 0.57272 28.60125 5.595 5.09e-06 ***
## condition1:group1 0.05207 0.60072 31.61980 0.087 0.931475
## condition1:exp1 0.80868 0.36940 38.61847 2.189 0.034695 *
## group1:exp1 2.84308 0.69995 37.01803 4.062 0.000243 ***
## accuracy1:group1 -0.25876 0.41212 821.92631 -0.628 0.530254
## accuracy1:exp1 0.28597 0.39196 974.52750 0.730 0.465816
## condition1:group1:exp1 1.42160 0.82438 26.74382 1.724 0.096167 .
## accuracy1:group1:exp1 1.23295 0.78442 976.38331 1.572 0.116320
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# no warnings
save(F02, file = "F02.Rdata")
#*Nice table output*---------------------------------
# sjPlot::tab_model((F02), show.se = TRUE, show.stat = TRUE, collapse.ci=TRUE, file="F02.doc")
# clean workspace ---------------------------------
# ls()
rm(corr5_F02, m0_F02, zcp_F02, zcp1_F02, zcp2_F02, zcp3_F02, zcp4_F02, zcp5_F02)4.4 Model reduction - pause (P)
# all "REML = FALSE" to use ANOVA model comparison later on
# INITIAL model with indicator variables ---------------------------------
summary(m0_P <- lmer(pause3 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd+cor.inc+
cor.inc_lhd.rhd|item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## pause3 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + cor.inc_rd.rpt + grpd.ungr_rd.rpt | Subj) +
## (1 + lhd.rhd + cor.inc + cor.inc_lhd.rhd | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8029.1 8297.5 -3962.6 7925.1 1237
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.5594 -0.3629 -0.0635 0.3633 6.6807
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 7.2997 2.7018
## cor.inc 4.3047 2.0748 -0.19
## grpd.ungr 8.5614 2.9260 0.22 0.52
## rd.rpt 12.8623 3.5864 0.27 -0.16 0.36
## grpd.ungr_cor.inc 50.7980 7.1273 0.78 0.07 0.09 -0.30
## cor.inc_rd.rpt 10.3345 3.2147 0.36 0.67 0.40 -0.46 0.76
## grpd.ungr_rd.rpt 28.3623 5.3256 0.30 0.10 0.74 0.69 0.05
## item (Intercept) 0.2771 0.5264
## lhd.rhd 0.1430 0.3782 0.08
## cor.inc 0.6806 0.8250 -0.97 0.15
## cor.inc_lhd.rhd 3.6581 1.9126 0.71 -0.64 -0.86
## Residual 22.7076 4.7652
##
##
##
##
##
##
##
## 0.07
##
##
##
##
##
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.9587 0.5407 31.8062 16.568 < 2e-16 ***
## grpd.ungr 8.2968 0.6928 31.1717 11.975 3.39e-13 ***
## lhd.rhd 2.0314 0.8598 31.2238 2.363 0.02455 *
## cor.inc 2.0552 0.6091 30.2512 3.374 0.00204 **
## rd.rpt 4.5842 0.8243 27.2811 5.561 6.54e-06 ***
## grpd.ungr_lhd.rhd 1.0749 1.3168 31.5545 0.816 0.42046
## grpd.ungr_cor.inc 17.7254 1.5717 31.6610 11.278 1.27e-12 ***
## cor.inc_lhd.rhd -0.2628 1.2298 28.8501 -0.214 0.83226
## grpd.ungr_rd.rpt 7.3703 1.2044 27.8586 6.120 1.36e-06 ***
## lhd.rhd_rd.rpt -0.7190 1.5242 28.8788 -0.472 0.64067
## cor.inc_rd.rpt 0.5114 1.1437 33.1015 0.447 0.65769
## grpd.ungr_lhd.rhd_rd.rpt -0.6326 2.3589 27.7217 -0.268 0.79054
## cor.inc_lhd.rhd_rd.rpt 2.0219 2.1337 36.4327 0.948 0.34956
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
# correlations seem okay for Subj and for items
summary(rePCA(m0_P))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5] [,6] [,7]
## Standard deviation 1.6691 1.3482 0.7525 0.38927 0.27237 0.0005812 0.0001294
## Proportion of Variance 0.5163 0.3369 0.1050 0.02808 0.01375 0.0000000 0.0000000
## Cumulative Proportion 0.5163 0.8532 0.9582 0.98625 1.00000 1.0000000 1.0000000
##
## $item
## Importance of components:
## [,1] [,2] [,3] [,4]
## Standard deviation 0.4398 0.12718 8.639e-05 0
## Proportion of Variance 0.9228 0.07718 0.000e+00 0
## Cumulative Proportion 0.9228 1.00000 1.000e+00 1# 5 of 7 components for Subj | 2 of 4 components for item contribute (variance)
# try model reduction of randEff ---------------------------------
## start with ZERO CORRELATION model =================================
summary(zcp_P <- lmer(pause3 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
cor.inc_rd.rpt+
grpd.ungr_rd.rpt||Subj)+
(1+lhd.rhd+cor.inc+
cor.inc_lhd.rhd||item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## pause3 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + cor.inc_rd.rpt + grpd.ungr_rd.rpt ||
## Subj) + (1 + lhd.rhd + cor.inc + cor.inc_lhd.rhd || item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8063.9 8192.9 -4006.9 8013.9 1264
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4266 -0.3456 -0.0304 0.3188 6.8826
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subj (Intercept) 7.849e+00 2.8015381
## Subj.1 cor.inc 4.978e+00 2.2311549
## Subj.2 grpd.ungr 9.762e+00 3.1244283
## Subj.3 rd.rpt 8.867e+00 2.9776854
## Subj.4 grpd.ungr_cor.inc 4.453e+01 6.6730943
## Subj.5 cor.inc_rd.rpt 7.430e-08 0.0002726
## Subj.6 grpd.ungr_rd.rpt 2.165e+01 4.6529683
## item (Intercept) 0.000e+00 0.0000000
## item.1 lhd.rhd 0.000e+00 0.0000000
## item.2 cor.inc 0.000e+00 0.0000000
## item.3 cor.inc_lhd.rhd 3.170e+00 1.7805591
## Residual 2.358e+01 4.8564295
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 9.0685 0.5564 31.7393 16.298 < 2e-16 ***
## grpd.ungr 8.5559 0.7138 34.9168 11.986 6.27e-14 ***
## lhd.rhd 2.5681 1.1139 31.8543 2.305 0.027808 *
## cor.inc 2.4926 0.6099 28.3101 4.087 0.000327 ***
## rd.rpt 4.3801 0.7172 32.7562 6.107 7.25e-07 ***
## grpd.ungr_lhd.rhd 1.7753 1.4111 34.0911 1.258 0.216909
## grpd.ungr_cor.inc 17.7410 1.4766 33.9227 12.015 8.99e-14 ***
## cor.inc_lhd.rhd 0.3076 1.2691 29.3077 0.242 0.810167
## grpd.ungr_rd.rpt 6.9718 1.1236 28.4198 6.205 9.91e-07 ***
## lhd.rhd_rd.rpt -1.0547 1.4312 32.4925 -0.737 0.466441
## cor.inc_rd.rpt 0.1846 0.8994 624.1574 0.205 0.837460
## grpd.ungr_lhd.rhd_rd.rpt -0.6249 2.2559 29.1007 -0.277 0.783745
## cor.inc_lhd.rhd_rd.rpt 1.3933 1.9449 88.5379 0.716 0.475638
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
'Random effects:
Groups Name Variance Std.Dev.
Subj (Intercept) 7.848 2.801
Subj.1 cor.inc 4.978 2.231
Subj.2 grpd.ungr 9.762 3.124
Subj.3 rd.rpt 8.867 2.978
Subj.4 grpd.ungr_cor.inc 44.529 6.673
Subj.5 cor.inc_rd.rpt 0.000 0.000 # no variance
Subj.6 grpd.ungr_rd.rpt 21.651 4.653
item (Intercept) 0.000 0.000
item.1 lhd.rhd 0.000 0.000 # no variance
item.2 cor.inc 0.000 0.000 # no variance
item.3 cor.inc_lhd.rhd 3.171 1.781
Residual 23.585 4.856
Number of obs: 1289, groups: Subj, 32; item, 24'## [1] "Random effects:\n Groups Name Variance Std.Dev.\n Subj (Intercept) 7.848 2.801 \n Subj.1 cor.inc 4.978 2.231 \n Subj.2 grpd.ungr 9.762 3.124 \n Subj.3 rd.rpt 8.867 2.978 \n Subj.4 grpd.ungr_cor.inc 44.529 6.673 \n Subj.5 cor.inc_rd.rpt 0.000 0.000 # no variance\n Subj.6 grpd.ungr_rd.rpt 21.651 4.653 \n item (Intercept) 0.000 0.000 \n item.1 lhd.rhd 0.000 0.000 # no variance \n item.2 cor.inc 0.000 0.000 # no variance\n item.3 cor.inc_lhd.rhd 3.171 1.781 \n Residual 23.585 4.856 \n Number of obs: 1289, groups: Subj, 32; item, 24"## ANOVA model comparison - initial model - zcp model =================================
anova(m0_P, zcp_P)# m0_P has signif. better fit - removing all correlations too harsh
## back to initial model =================================
# try exclusion of randEff slopes with regards to zcp
### FIRST #############################
# try exclusion of "cor.inc_rd.rpt" as slope per Subj
summary(m1_P <- lmer(pause3 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd+cor.inc+
cor.inc_lhd.rhd|item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## pause3 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + grpd.ungr_rd.rpt | Subj) + (1 + lhd.rhd +
## cor.inc + cor.inc_lhd.rhd | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8024.4 8256.7 -3967.2 7934.4 1244
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.5182 -0.3580 -0.0663 0.3493 6.6451
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 7.80079 2.7930
## cor.inc 5.09424 2.2570 -0.27
## grpd.ungr 9.83907 3.1367 0.25 0.61
## rd.rpt 11.06323 3.3261 0.39 0.01 0.44
## grpd.ungr_cor.inc 47.12526 6.8648 0.82 -0.11 0.15 -0.12
## grpd.ungr_rd.rpt 25.51835 5.0516 0.33 -0.01 0.65 0.80 0.02
## item (Intercept) 0.19577 0.4425
## lhd.rhd 0.02583 0.1607 -0.41
## cor.inc 0.58776 0.7667 -0.99 0.52
## cor.inc_lhd.rhd 3.08416 1.7562 0.96 -0.66 -0.99
## Residual 23.02807 4.7988
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.949921 0.548281 31.719094 16.324 < 2e-16 ***
## grpd.ungr 8.318558 0.713032 33.182192 11.666 2.76e-13 ***
## lhd.rhd 1.949930 0.837642 33.498895 2.328 0.0261 *
## cor.inc 2.002837 0.626261 27.133272 3.198 0.0035 **
## rd.rpt 4.823518 0.764060 33.070591 6.313 3.82e-07 ***
## grpd.ungr_lhd.rhd 1.375990 1.380263 32.592761 0.997 0.3262
## grpd.ungr_cor.inc 17.611972 1.527892 33.443383 11.527 3.42e-13 ***
## cor.inc_lhd.rhd 0.157974 1.266159 25.962715 0.125 0.9017
## grpd.ungr_rd.rpt 7.327576 1.169706 27.757274 6.264 9.36e-07 ***
## lhd.rhd_rd.rpt -0.785435 1.459946 33.487118 -0.538 0.5941
## cor.inc_rd.rpt -0.007487 0.958254 89.913359 -0.008 0.9938
## grpd.ungr_lhd.rhd_rd.rpt -0.486638 2.292581 28.667437 -0.212 0.8334
## cor.inc_lhd.rhd_rd.rpt 2.563852 1.946562 68.057792 1.317 0.1922
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
summary(rePCA(m1_P))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5] [,6]
## Standard deviation 1.5234 1.2950 0.65188 0.39631 0.20588 0
## Proportion of Variance 0.5021 0.3628 0.09194 0.03398 0.00917 0
## Cumulative Proportion 0.5021 0.8649 0.95685 0.99083 1.00000 1
##
## $item
## Importance of components:
## [,1] [,2] [,3] [,4]
## Standard deviation 0.4090 0.04284 8.107e-05 5.619e-06
## Proportion of Variance 0.9891 0.01085 0.000e+00 0.000e+00
## Cumulative Proportion 0.9891 1.00000 1.000e+00 1.000e+00# 1 of 6 components for Subj | 2 of 4 components for item have zero variance
### ANOVA model comparison - initial model - reduced model #############################
anova(m0_P, m1_P)# not significant and AIC smaller for m1_P - preferred
### SECOND #############################
# try exclusion also of "lhd.rhd" as slope per item
summary(m2_P <- lmer(pause3 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
grpd.ungr_rd.rpt|Subj)+
(1+cor.inc+
cor.inc_lhd.rhd|item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## pause3 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + grpd.ungr_rd.rpt | Subj) + (1 + cor.inc +
## cor.inc_lhd.rhd | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8016.4 8228.1 -3967.2 7934.4 1248
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.5196 -0.3589 -0.0665 0.3494 6.6648
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 7.8039 2.7935
## cor.inc 5.0745 2.2527 -0.27
## grpd.ungr 9.8340 3.1359 0.25 0.61
## rd.rpt 11.0756 3.3280 0.39 0.01 0.44
## grpd.ungr_cor.inc 47.0909 6.8623 0.82 -0.11 0.15 -0.12
## grpd.ungr_rd.rpt 25.5208 5.0518 0.33 0.00 0.65 0.80 0.02
## item (Intercept) 0.2083 0.4565
## cor.inc 0.6231 0.7893 -1.00
## cor.inc_lhd.rhd 2.6435 1.6259 1.00 -1.00
## Residual 23.0385 4.7998
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.948821 0.548837 31.952675 16.305 < 2e-16 ***
## grpd.ungr 8.310981 0.714354 35.258137 11.634 1.26e-13 ***
## lhd.rhd 1.942741 0.834519 34.799515 2.328 0.02585 *
## cor.inc 2.001122 0.626867 27.475175 3.192 0.00352 **
## rd.rpt 4.820944 0.765598 34.588301 6.297 3.30e-07 ***
## grpd.ungr_lhd.rhd 1.352954 1.378222 33.695812 0.982 0.33326
## grpd.ungr_cor.inc 17.619918 1.528850 33.496345 11.525 3.36e-13 ***
## cor.inc_lhd.rhd 0.149601 1.257615 27.813327 0.119 0.90617
## grpd.ungr_rd.rpt 7.344029 1.170218 28.759715 6.276 7.77e-07 ***
## lhd.rhd_rd.rpt -0.775643 1.458619 34.119850 -0.532 0.59833
## cor.inc_rd.rpt 0.005364 0.961115 93.109017 0.006 0.99556
## grpd.ungr_lhd.rhd_rd.rpt -0.417663 2.294091 29.058527 -0.182 0.85680
## cor.inc_lhd.rhd_rd.rpt 2.573566 1.926683 123.584475 1.336 0.18409
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
# Model failed to converge with 1 negative eigenvalue: -1.7e+01
### THIRD #############################
# back to FIRST REDUCED model
# try exclusion also of "cor.inc" as slope per item
summary(m3_P <- lmer(pause3 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd+
cor.inc_lhd.rhd|item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## pause3 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + grpd.ungr_rd.rpt | Subj) + (1 + lhd.rhd +
## cor.inc_lhd.rhd | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8019.0 8230.6 -3968.5 7937.0 1248
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.5110 -0.3505 -0.0638 0.3542 6.8294
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 7.76116 2.7859
## cor.inc 5.37409 2.3182 -0.25
## grpd.ungr 9.96203 3.1563 0.26 0.59
## rd.rpt 10.84624 3.2934 0.40 0.01 0.45
## grpd.ungr_cor.inc 47.23947 6.8731 0.82 -0.09 0.16 -0.11
## grpd.ungr_rd.rpt 25.32193 5.0321 0.35 -0.05 0.65 0.80 0.02
## item (Intercept) 0.02525 0.1589
## lhd.rhd 0.10694 0.3270 -0.84
## cor.inc_lhd.rhd 4.65239 2.1569 0.94 -0.97
## Residual 23.17005 4.8135
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.9642 0.5409 30.6633 16.574 < 2e-16 ***
## grpd.ungr 8.4258 0.6961 33.1977 12.104 1.01e-13 ***
## lhd.rhd 1.9457 0.8511 34.9442 2.286 0.02843 *
## cor.inc 2.0548 0.6152 25.0731 3.340 0.00262 **
## rd.rpt 4.7831 0.7410 33.6425 6.455 2.33e-07 ***
## grpd.ungr_lhd.rhd 1.4307 1.3963 33.2587 1.025 0.31291
## grpd.ungr_cor.inc 17.5383 1.4858 31.2194 11.804 4.80e-13 ***
## cor.inc_lhd.rhd 0.2253 1.3092 26.0353 0.172 0.86470
## grpd.ungr_rd.rpt 7.0596 1.1576 27.6135 6.099 1.49e-06 ***
## lhd.rhd_rd.rpt -0.7651 1.4584 34.4065 -0.525 0.60321
## cor.inc_rd.rpt -0.1191 0.9099 345.4753 -0.131 0.89598
## grpd.ungr_lhd.rhd_rd.rpt -0.6383 2.2883 28.8280 -0.279 0.78227
## cor.inc_lhd.rhd_rd.rpt 2.5114 2.0241 53.1778 1.241 0.22014
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
# Model failed to converge with max|grad| = 0.0207696 (tol = 0.002, component 1)
### FOURTH #############################
# back to FIRST REDUCED model
# try exclusion also of "cor.inc_lhd.rhd" as slope per item
summary(m4_P <- lmer(pause3 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd+cor.inc|item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## pause3 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + grpd.ungr_rd.rpt | Subj) + (1 + lhd.rhd +
## cor.inc | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8019.5 8231.1 -3968.7 7937.5 1248
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4941 -0.3618 -0.0633 0.3446 6.8902
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 7.6740 2.7702
## cor.inc 5.0733 2.2524 -0.27
## grpd.ungr 9.8985 3.1462 0.24 0.62
## rd.rpt 11.2686 3.3569 0.39 0.01 0.43
## grpd.ungr_cor.inc 46.9611 6.8528 0.82 -0.11 0.15 -0.12
## grpd.ungr_rd.rpt 25.2664 5.0266 0.34 0.01 0.65 0.80 0.03
## item (Intercept) 0.3274 0.5722
## lhd.rhd 0.2040 0.4517 1.00
## cor.inc 0.9180 0.9581 -1.00 -1.00
## Residual 23.1687 4.8134
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.93784 0.55072 32.71151 16.229 < 2e-16 ***
## grpd.ungr 8.28697 0.73114 34.69898 11.334 3.26e-13 ***
## lhd.rhd 1.87398 0.82500 33.88069 2.271 0.02959 *
## cor.inc 2.03715 0.63883 27.47368 3.189 0.00355 **
## rd.rpt 4.77708 0.78375 34.26532 6.095 6.29e-07 ***
## grpd.ungr_lhd.rhd 1.21257 1.39368 34.02248 0.870 0.39037
## grpd.ungr_cor.inc 17.68798 1.53287 32.27654 11.539 5.43e-13 ***
## cor.inc_lhd.rhd 0.12183 1.21667 26.27744 0.100 0.92100
## grpd.ungr_rd.rpt 7.32037 1.17762 28.78530 6.216 9.10e-07 ***
## lhd.rhd_rd.rpt -0.68361 1.47888 34.35887 -0.462 0.64682
## cor.inc_rd.rpt 0.08517 0.99386 58.35029 0.086 0.93200
## grpd.ungr_lhd.rhd_rd.rpt -0.27075 2.29671 29.12705 -0.118 0.90697
## cor.inc_lhd.rhd_rd.rpt 2.66012 1.81970 355.05148 1.462 0.14467
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
summary(rePCA(m4_P))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5] [,6]
## Standard deviation 1.5147 1.2897 0.65459 0.38981 0.20764 9.061e-06
## Proportion of Variance 0.5008 0.3631 0.09353 0.03317 0.00941 0.000e+00
## Cumulative Proportion 0.5008 0.8639 0.95742 0.99059 1.00000 1.000e+00
##
## $item
## Importance of components:
## [,1] [,2] [,3]
## Standard deviation 0.2501 1.335e-06 0
## Proportion of Variance 1.0000 0.000e+00 0
## Cumulative Proportion 1.0000 1.000e+00 1# 1 of 6 components for Subj | 2 of 3 components for item have zero variance
### ANOVA model comparison - initial model - reduced model #############################
anova(m0_P, m4_P)# not significant and AIC smaller for m4_P - preferred
### FIFTH #############################
# try exclusion also of "lhd.rhd" as slope per item
summary(m5_P <- lmer(pause3 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
grpd.ungr_rd.rpt|Subj)+
(1+cor.inc|item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## pause3 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + grpd.ungr_rd.rpt | Subj) + (1 + cor.inc | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8015.0 8211.1 -3969.5 7939.0 1251
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4806 -0.3661 -0.0631 0.3485 6.8999
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 7.6499 2.7658
## cor.inc 5.0889 2.2558 -0.27
## grpd.ungr 9.8426 3.1373 0.23 0.63
## rd.rpt 11.1125 3.3335 0.39 0.01 0.44
## grpd.ungr_cor.inc 47.3394 6.8804 0.82 -0.11 0.15 -0.12
## grpd.ungr_rd.rpt 24.9721 4.9972 0.34 0.00 0.65 0.80 0.03
## item (Intercept) 0.1878 0.4333
## cor.inc 0.7248 0.8514 -1.00
## Residual 23.2779 4.8247
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.93014 0.54435 31.68665 16.405 < 2e-16 ***
## grpd.ungr 8.29449 0.71318 33.07472 11.630 3.14e-13 ***
## lhd.rhd 1.86581 0.81542 33.58536 2.288 0.02855 *
## cor.inc 2.07624 0.63212 26.26237 3.285 0.00289 **
## rd.rpt 4.75769 0.76472 32.81165 6.221 5.16e-07 ***
## grpd.ungr_lhd.rhd 1.24228 1.37846 33.47333 0.901 0.37391
## grpd.ungr_cor.inc 17.72467 1.52533 31.54350 11.620 6.22e-13 ***
## cor.inc_lhd.rhd 0.15660 1.21608 25.91478 0.129 0.89853
## grpd.ungr_rd.rpt 7.17001 1.15699 27.94072 6.197 1.09e-06 ***
## lhd.rhd_rd.rpt -0.71160 1.46007 33.94266 -0.487 0.62913
## cor.inc_rd.rpt 0.07505 0.97564 57.36875 0.077 0.93895
## grpd.ungr_lhd.rhd_rd.rpt -0.64562 2.27537 28.97303 -0.284 0.77862
## cor.inc_lhd.rhd_rd.rpt 2.67822 1.81616 349.27549 1.475 0.14120
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
summary(rePCA(m5_P))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5] [,6]
## Standard deviation 1.5171 1.2790 0.65337 0.38784 0.19747 3.341e-06
## Proportion of Variance 0.5054 0.3592 0.09374 0.03303 0.00856 0.000e+00
## Cumulative Proportion 0.5054 0.8647 0.95841 0.99144 1.00000 1.000e+00
##
## $item
## Importance of components:
## [,1] [,2]
## Standard deviation 0.198 7.709e-06
## Proportion of Variance 1.000 0.000e+00
## Cumulative Proportion 1.000 1.000e+00# 1 of 6 components for Subj | 1 of 2 components for item have zero variance
### ANOVA model comparison - m4_P - m5_P #############################
anova(m4_P, m5_P)# not significant and AIC smaller for m5_P - preferred
### SIXTH #############################
# try exclusion also of "cor.inc" as slope per item - intercept only
summary(m6_P <- lmer(pause3 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
grpd.ungr_rd.rpt|Subj)+
(1|item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## pause3 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + grpd.ungr_rd.rpt | Subj) + (1 | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8011.9 8197.8 -3970.0 7939.9 1253
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4598 -0.3677 -0.0655 0.3424 7.0755
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 7.702e+00 2.77532
## cor.inc 5.089e+00 2.25584 -0.27
## grpd.ungr 9.979e+00 3.15892 0.25 0.61
## rd.rpt 1.094e+01 3.30787 0.41 0.02 0.46
## grpd.ungr_cor.inc 4.674e+01 6.83637 0.82 -0.12 0.15 -0.11
## grpd.ungr_rd.rpt 2.493e+01 4.99277 0.35 -0.01 0.66 0.81 0.03
## item (Intercept) 1.211e-08 0.00011
## Residual 2.347e+01 4.84438
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.95019 0.53759 30.70565 16.649 < 2e-16 ***
## grpd.ungr 8.38150 0.69131 33.58594 12.124 8.13e-14 ***
## lhd.rhd 1.89765 0.81713 34.17574 2.322 0.02630 *
## cor.inc 2.07156 0.60665 25.81942 3.415 0.00212 **
## rd.rpt 4.74578 0.73918 33.92124 6.420 2.49e-07 ***
## grpd.ungr_lhd.rhd 1.31088 1.38334 33.55847 0.948 0.35010
## grpd.ungr_cor.inc 17.62183 1.47667 31.20921 11.933 3.64e-13 ***
## cor.inc_lhd.rhd 0.16200 1.21488 25.91339 0.133 0.89495
## grpd.ungr_rd.rpt 6.95256 1.14014 28.10458 6.098 1.39e-06 ***
## lhd.rhd_rd.rpt -0.65959 1.45254 34.25120 -0.454 0.65262
## cor.inc_rd.rpt -0.03247 0.90747 351.22544 -0.036 0.97148
## grpd.ungr_lhd.rhd_rd.rpt -0.53740 2.27302 29.06302 -0.236 0.81476
## cor.inc_lhd.rhd_rd.rpt 2.56241 1.81035 353.10480 1.415 0.15783
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
# Model failed to converge with 1 negative eigenvalue: -6.1e+00
### SEVENTH #############################
# BACK to m4_P
# try exclusion of "cor.inc" as slope per item instead of "lhd.rhd"
summary(m7_P <- lmer(pause3 ~
1+grpd.ungr+lhd.rhd+cor.inc+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_cor.inc+
cor.inc_lhd.rhd+grpd.ungr_rd.rpt+
lhd.rhd_rd.rpt+cor.inc_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+cor.inc_lhd.rhd_rd.rpt+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd|item),
REML = FALSE,pDat_mm))## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula:
## pause3 ~ 1 + grpd.ungr + lhd.rhd + cor.inc + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_cor.inc + cor.inc_lhd.rhd + grpd.ungr_rd.rpt +
## lhd.rhd_rd.rpt + cor.inc_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## cor.inc_lhd.rhd_rd.rpt + (1 + cor.inc + grpd.ungr + rd.rpt +
## grpd.ungr_cor.inc + grpd.ungr_rd.rpt | Subj) + (1 + lhd.rhd | item)
## Data: pDat_mm
##
## AIC BIC logLik deviance df.resid
## 8015.7 8211.9 -3969.9 7939.7 1251
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4693 -0.3496 -0.0644 0.3397 7.0799
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 7.71883 2.7783
## cor.inc 5.07372 2.2525 -0.27
## grpd.ungr 10.02589 3.1664 0.25 0.61
## rd.rpt 10.99376 3.3157 0.40 0.01 0.45
## grpd.ungr_cor.inc 46.57779 6.8248 0.82 -0.12 0.15 -0.11
## grpd.ungr_rd.rpt 25.07541 5.0075 0.35 -0.02 0.66 0.81 0.03
## item (Intercept) 0.04957 0.2226
## lhd.rhd 0.08933 0.2989 1.00
## Residual 23.40518 4.8379
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.95525 0.54003 30.91276 16.583 < 2e-16 ***
## grpd.ungr 8.38203 0.69838 33.25103 12.002 1.24e-13 ***
## lhd.rhd 1.90989 0.82162 34.10226 2.325 0.02619 *
## cor.inc 2.04901 0.60650 25.92744 3.378 0.00231 **
## rd.rpt 4.75765 0.74604 33.68425 6.377 2.92e-07 ***
## grpd.ungr_lhd.rhd 1.30691 1.39098 33.54422 0.940 0.35416
## grpd.ungr_cor.inc 17.59049 1.47540 31.22239 11.923 3.71e-13 ***
## cor.inc_lhd.rhd 0.13121 1.21458 26.04676 0.108 0.91480
## grpd.ungr_rd.rpt 6.97248 1.15675 27.47555 6.028 1.83e-06 ***
## lhd.rhd_rd.rpt -0.64188 1.45993 34.24185 -0.440 0.66294
## cor.inc_rd.rpt -0.04458 0.90792 353.73071 -0.049 0.96086
## grpd.ungr_lhd.rhd_rd.rpt -0.51727 2.29016 29.03852 -0.226 0.82289
## cor.inc_lhd.rhd_rd.rpt 2.52673 1.81156 355.09139 1.395 0.16395
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
summary(rePCA(m7_P))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4] [,5] [,6]
## Standard deviation 1.5044 1.2790 0.64670 0.38973 0.19191 0
## Proportion of Variance 0.5023 0.3630 0.09281 0.03371 0.00817 0
## Cumulative Proportion 0.5023 0.8653 0.95812 0.99183 1.00000 1
##
## $item
## Importance of components:
## [,1] [,2]
## Standard deviation 0.07704 1.904e-05
## Proportion of Variance 1.00000 0.000e+00
## Cumulative Proportion 1.00000 1.000e+00# 1 of 6 components for Subj | 1 of 2 components for item have zero variance
# same for m5_P
### ANOVA model comparison - initial model - m7_P #############################
anova(m0_P, m7_P)# not significant and AIC smaller for m7_P - preferred
### ANOVA model comparison - initial model - m5_P #############################
anova(m0_P, m5_P)# not significant and AIC smaller for m7_P - preferred
# both ANOVAs obtain almost identical measures
# stay with m7_P
# FINAL model ---------------------------------
## model with factors and * for fixed effects =================================
# use REML=TRUE
summary(P <- lmer(pause3 ~
condition*accuracy+condition*group*exp+
group*exp*accuracy+
(1+cor.inc+grpd.ungr+rd.rpt+
grpd.ungr_cor.inc+
grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd|item),
REML = TRUE,pDat_mm))## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: pause3 ~ condition * accuracy + condition * group * exp + group *
## exp * accuracy + (1 + cor.inc + grpd.ungr + rd.rpt + grpd.ungr_cor.inc +
## grpd.ungr_rd.rpt | Subj) + (1 + lhd.rhd | item)
## Data: pDat_mm
##
## REML criterion at convergence: 7918.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.4809 -0.3452 -0.0567 0.3338 7.0880
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 8.25955 2.8739
## cor.inc 6.02594 2.4548 -0.27
## grpd.ungr 11.16365 3.3412 0.26 0.57
## rd.rpt 12.10076 3.4786 0.40 0.01 0.44
## grpd.ungr_cor.inc 48.73116 6.9808 0.81 -0.10 0.14 -0.12
## grpd.ungr_rd.rpt 28.19503 5.3099 0.34 -0.06 0.63 0.80 0.01
## item (Intercept) 0.08123 0.2850
## lhd.rhd 0.14436 0.3799 1.00
## Residual 23.41247 4.8386
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.9745 0.5588 29.7961 16.060 3.25e-16 ***
## condition1 8.4361 0.7327 31.6435 11.514 7.54e-13 ***
## accuracy1 2.0183 0.6377 23.9710 3.165 0.00418 **
## group1 1.9202 0.8571 31.7202 2.240 0.03220 *
## exp1 4.8107 0.7782 31.9286 6.182 6.48e-07 ***
## condition1:accuracy1 17.4744 1.5081 29.9484 11.587 1.36e-12 ***
## condition1:group1 1.3088 1.4557 31.5896 0.899 0.37541
## condition1:exp1 6.9831 1.2164 26.3545 5.741 4.58e-06 ***
## group1:exp1 -0.5902 1.5188 31.9723 -0.389 0.70016
## accuracy1:group1 0.1183 1.2768 24.2112 0.093 0.92696
## accuracy1:exp1 -0.1541 0.9227 362.3933 -0.167 0.86743
## condition1:group1:exp1 -0.4799 2.3977 27.5151 -0.200 0.84282
## accuracy1:group1:exp1 2.4604 1.8416 362.2902 1.336 0.18238
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
save(P, file = "P.Rdata")
# *Nice table output* ---------------------------------
# sjPlot::tab_model((P), show.se = TRUE, show.stat = TRUE, collapse.ci=TRUE, file="P.doc")
# clean workspace ---------------------------------
# ls()
rm(m0_P, m1_P, m2_P, m3_P, m4_P, m5_P, m6_P, m7_P, zcp_P)5 Final lengthening model (FL2) for name2
5.1 Visualizations and Post-Hoc-Tests
# interaction visualization ---------------------------------
# condition:accuracy =================================
### condition ~ accuracy #############################
load("FL2.Rdata")
emmeans::emmip(FL2, condition ~ accuracy, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - % final lengthening on name2 \n(relative to duration of name2)",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
### accuracy ~ condition #############################
emmeans::emmip(FL2, accuracy ~ condition, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - % final lengthening on name2 \n(relative to duration of name2)",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
# POST-HOC testing ---------------------------------
# what are the significant differences
emmeans::emmeans(FL2, pairwise ~ condition*accuracy) %>%
summary(by = NULL, adjust = "bonferroni")## $emmeans
## condition accuracy emmean SE df lower.CL upper.CL
## grouped correct 48.1 1.21 35.7 44.9 51.3
## ungrouped correct 38.7 1.48 49.6 34.8 42.5
## grouped incorrect 43.5 1.19 36.1 40.4 46.6
## ungrouped incorrect 43.5 1.40 37.4 39.9 47.2
##
## Results are averaged over the levels of: group, exp
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 4 estimates
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## grouped correct - ungrouped correct 9.4548 1.69 39.3 5.607 <.0001
## grouped correct - grouped incorrect 4.6467 1.02 29.4 4.547 0.0005
## grouped correct - ungrouped incorrect 4.5924 1.53 42.1 2.993 0.0276
## ungrouped correct - grouped incorrect -4.8081 1.49 38.0 -3.227 0.0155
## ungrouped correct - ungrouped incorrect -4.8624 1.23 34.1 -3.948 0.0022
## grouped incorrect - ungrouped incorrect -0.0543 1.46 32.0 -0.037 1.0000
##
## Results are averaged over the levels of: group, exp
## Degrees-of-freedom method: kenward-roger
## P value adjustment: bonferroni method for 6 tests# interaction visualization ---------------------------------
# group:exp =================================
### group ~ exp #############################
emmeans::emmip(FL2, group ~ exp, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - % final lengthening on name2 \n(relative to duration of name2)",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
### exp ~ group #############################
emmeans::emmip(FL2, exp ~ group, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - % final lengthening on name2 \n(relative to duration of name2)",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
# POST-HOC testing ---------------------------------
# what are the significant differences
emmeans::emmeans(FL2, pairwise ~ group*exp) %>%
summary(by = NULL, adjust = "bonferroni")## $emmeans
## group exp emmean SE df lower.CL upper.CL
## LHDP reading_aloud 45.3 1.49 37.7 41.3 49.2
## RHDP reading_aloud 43.4 1.67 44.2 39.1 47.7
## LHDP repetition 45.7 1.54 34.7 41.6 49.7
## RHDP repetition 39.6 1.65 42.3 35.2 43.9
##
## Results are averaged over the levels of: condition, accuracy
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 4 estimates
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## LHDP reading_aloud - RHDP reading_aloud 1.852 1.99 37.1 0.932 1.0000
## LHDP reading_aloud - LHDP repetition -0.405 1.51 38.3 -0.268 1.0000
## LHDP reading_aloud - RHDP repetition 5.700 2.25 58.6 2.538 0.0829
## RHDP reading_aloud - LHDP repetition -2.257 2.29 59.4 -0.987 1.0000
## RHDP reading_aloud - RHDP repetition 3.848 1.90 28.6 2.021 0.3163
## LHDP repetition - RHDP repetition 6.104 2.03 34.6 3.013 0.0289
##
## Results are averaged over the levels of: condition, accuracy
## Degrees-of-freedom method: kenward-roger
## P value adjustment: bonferroni method for 6 tests6 F0_range model (F02) for name2
6.1 Visualizations and Post-Hoc-Tests
# interaction visualization ---------------------------------
load("F02.Rdata")
# condition:accuracy =================================
### condition ~ accuracy #############################
emmeans::emmip(F02, condition ~ accuracy,
CIs = TRUE, CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - f0 range on name2 in semitones",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
### accuracy ~ condition #############################
emmeans::emmip(F02, accuracy ~ condition,
CIs = TRUE, CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - f0 range on name2 in semitones",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
# POST-HOC testing ---------------------------------
# what are the significant differences
emmeans::emmeans(F02, pairwise ~ condition*accuracy) %>%
summary(by = NULL, adjust = "bonferroni")## $emmeans
## condition accuracy emmean SE df lower.CL upper.CL
## grouped correct 7.29 0.500 33.1 5.97 8.61
## ungrouped correct 3.62 0.289 32.9 2.86 4.38
## grouped incorrect 5.30 0.452 38.8 4.12 6.49
## ungrouped incorrect 4.84 0.451 39.4 3.66 6.02
##
## Results are averaged over the levels of: group, exp
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 4 estimates
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## grouped correct - ungrouped correct 3.666 0.374 33.7 9.790 <.0001
## grouped correct - grouped incorrect 1.983 0.380 71.7 5.224 <.0001
## grouped correct - ungrouped incorrect 2.444 0.370 75.2 6.601 <.0001
## ungrouped correct - grouped incorrect -1.683 0.355 51.6 -4.743 0.0001
## ungrouped correct - ungrouped incorrect -1.222 0.343 58.6 -3.562 0.0044
## grouped incorrect - ungrouped incorrect 0.461 0.449 28.6 1.028 1.0000
##
## Results are averaged over the levels of: group, exp
## Degrees-of-freedom method: kenward-roger
## P value adjustment: bonferroni method for 6 tests# interaction visualization ---------------------------------
# condition:exp =================================
### condition ~ exp #############################
emmeans::emmip(F02, condition ~ exp,
CIs = TRUE, CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - f0 range on name2 in semitones",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
### exp ~ condition #############################
emmeans::emmip(F02, exp ~ condition,
CIs = TRUE, CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - f0 range on name2 in semitones",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
# POST-HOC testing ---------------------------------
# what are the significant differences
emmeans::emmeans(F02, pairwise ~ condition*exp) %>%
summary(by = NULL, adjust = "bonferroni")## $emmeans
## condition exp emmean SE df lower.CL upper.CL
## grouped reading_aloud 7.69 0.479 36.6 6.43 8.95
## ungrouped reading_aloud 5.22 0.389 38.7 4.20 6.24
## grouped repetition 4.91 0.481 33.1 3.63 6.18
## ungrouped repetition 3.25 0.390 34.2 2.22 4.27
##
## Results are averaged over the levels of: accuracy, group
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 4 estimates
##
## $contrasts
## contrast estimate SE df t.ratio
## grouped reading_aloud - ungrouped reading_aloud 2.468 0.341 36.6 7.227
## grouped reading_aloud - grouped repetition 2.781 0.396 43.5 7.016
## grouped reading_aloud - ungrouped repetition 4.441 0.451 33.4 9.854
## ungrouped reading_aloud - grouped repetition 0.313 0.454 32.0 0.690
## ungrouped reading_aloud - ungrouped repetition 1.973 0.388 42.4 5.082
## grouped repetition - ungrouped repetition 1.659 0.353 40.4 4.696
## p.value
## <.0001
## <.0001
## <.0001
## 1.0000
## <.0001
## 0.0002
##
## Results are averaged over the levels of: accuracy, group
## Degrees-of-freedom method: kenward-roger
## P value adjustment: bonferroni method for 6 tests# interaction visualization ---------------------------------
# group:exp =================================
### exp ~ group #############################
emmeans::emmip(F02, exp ~ group, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - f0 range on name2 in semitones",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
### group ~ exp #############################
emmeans::emmip(F02, group ~ exp, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - f0 range on name2 in semitones",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
# POST-HOC testing ---------------------------------
# what are the significant differences
emmeans::emmeans(F02, pairwise ~ group*exp) %>%
summary(by = NULL, adjust = "bonferroni")## $emmeans
## group exp emmean SE df lower.CL upper.CL
## LHDP reading_aloud 7.63 0.597 38.7 6.07 9.20
## RHDP reading_aloud 5.27 0.535 39.2 3.87 6.67
## LHDP repetition 3.83 0.593 32.3 2.26 5.40
## RHDP repetition 4.32 0.519 34.8 2.95 5.69
##
## Results are averaged over the levels of: condition, accuracy
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 4 estimates
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## LHDP reading_aloud - RHDP reading_aloud 2.356 0.800 39.1 2.944 0.0326
## LHDP reading_aloud - LHDP repetition 3.799 0.532 34.3 7.140 <.0001
## LHDP reading_aloud - RHDP repetition 3.312 0.778 51.6 4.258 0.0005
## RHDP reading_aloud - LHDP repetition 1.443 0.779 47.8 1.851 0.4218
## RHDP reading_aloud - RHDP repetition 0.956 0.470 34.3 2.032 0.2998
## LHDP repetition - RHDP repetition -0.487 0.775 31.6 -0.629 1.0000
##
## Results are averaged over the levels of: condition, accuracy
## Degrees-of-freedom method: kenward-roger
## P value adjustment: bonferroni method for 6 tests7 Pause model (P) for pause after name2
7.1 Visualizations and Post-Hoc-Tests
# interaction visualization ---------------------------------
load("P.Rdata")
# condition:accuracy =================================
### condition ~ accuracy #############################
emmeans::emmip(P, condition ~ accuracy,
CIs = TRUE, CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - % pause after name2 \n(relative to utterance duration)",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
### accuracy ~ condition #############################
emmeans::emmip(P, accuracy ~ condition, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - % pause after name2 \n(relative to utterance duration)",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
# POST-HOC testing ---------------------------------
# what are the significant differences
emmeans::emmeans(P, pairwise ~ condition * accuracy) %>%
summary(by = NULL, adjust = "bonferroni")## $emmeans
## condition accuracy emmean SE df lower.CL upper.CL
## grouped correct 18.57 1.036 29.4 15.81 21.33
## ungrouped correct 1.40 0.358 27.0 0.44 2.35
## grouped incorrect 7.81 0.774 31.1 5.76 9.87
## ungrouped incorrect 8.12 1.056 28.9 5.30 10.93
##
## Results are averaged over the levels of: group, exp
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 4 estimates
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## grouped correct - ungrouped correct 17.173 1.055 29.6 16.278 <.0001
## grouped correct - grouped incorrect 10.755 1.053 28.8 10.217 <.0001
## grouped correct - ungrouped incorrect 10.454 1.183 28.5 8.835 <.0001
## ungrouped correct - grouped incorrect -6.418 0.802 31.2 -7.999 <.0001
## ungrouped correct - ungrouped incorrect -6.719 0.999 28.5 -6.726 <.0001
## grouped incorrect - ungrouped incorrect -0.301 1.116 30.3 -0.270 1.0000
##
## Results are averaged over the levels of: group, exp
## Degrees-of-freedom method: kenward-roger
## P value adjustment: bonferroni method for 6 tests# interaction visualization ---------------------------------
# condition:exp =================================
### condition ~ exp #############################
emmeans::emmip(P, condition ~ exp, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - % pause after name2 \n(relative to utterance duration)",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
### exp ~ condition #############################
emmeans::emmip(P, exp ~ condition, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - % pause after name2 \n(relative to utterance duration)",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
# POST-HOC testing ---------------------------------
# what are the significant differences
emmeans::emmeans(P, pairwise ~ condition * exp) %>%
summary(by = NULL, adjust = "bonferroni")## $emmeans
## condition exp emmean SE df lower.CL upper.CL
## grouped reading_aloud 17.34 1.165 31.7 14.26 20.43
## ungrouped reading_aloud 5.42 0.729 32.3 3.49 7.34
## grouped repetition 9.04 0.744 27.3 7.05 11.03
## ungrouped repetition 4.10 0.674 30.7 2.31 5.88
##
## Results are averaged over the levels of: accuracy, group
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 4 estimates
##
## $contrasts
## contrast estimate SE df t.ratio
## grouped reading_aloud - ungrouped reading_aloud 11.93 1.131 29.6 10.547
## grouped reading_aloud - grouped repetition 8.30 1.258 29.4 6.601
## grouped reading_aloud - ungrouped repetition 13.25 1.275 35.4 10.390
## ungrouped reading_aloud - grouped repetition -3.63 0.897 29.7 -4.041
## ungrouped reading_aloud - ungrouped repetition 1.32 0.695 38.0 1.899
## grouped repetition - ungrouped repetition 4.94 0.807 26.2 6.124
## p.value
## <.0001
## <.0001
## <.0001
## 0.0021
## 0.3911
## <.0001
##
## Results are averaged over the levels of: accuracy, group
## Degrees-of-freedom method: kenward-roger
## P value adjustment: bonferroni method for 6 tests8 Model reduction - production accuracy
# INITIAL MAX model ---------------------------------
summary(rateM_i <- glmer(rate ~ condition*group*exp+
(1+condition*exp|Subj)+
(1+group|item),
family="binomial",
control=glmerControl(optCtrl=list(maxfun=1e5),
optimizer="bobyqa"),
pDat))## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: rate ~ condition * group * exp + (1 + condition * exp | Subj) +
## (1 + group | item)
## Data: pDat
## Control: glmerControl(optCtrl = list(maxfun = 1e+05), optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 1172.7 1281.1 -565.4 1130.7 1268
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5739 -0.1593 0.2867 0.4828 3.0638
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 2.57116 1.6035
## condition1 8.23428 2.8695 0.95
## exp1 4.13984 2.0347 0.34 0.38
## condition1:exp1 4.38362 2.0937 0.11 0.23 0.95
## item (Intercept) 0.01638 0.1280
## group1 0.13368 0.3656 -1.00
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.5818 0.3629 4.359 1.31e-05 ***
## condition1 0.1798 0.6630 0.271 0.786291
## group1 0.7818 0.6713 1.165 0.244157
## exp1 2.3320 0.5704 4.089 4.34e-05 ***
## condition1:group1 1.1029 1.2059 0.915 0.360408
## condition1:exp1 3.5188 0.9340 3.767 0.000165 ***
## group1:exp1 1.4909 1.0070 1.481 0.138704
## condition1:group1:exp1 2.2744 1.5297 1.487 0.137071
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) cndtn1 group1 exp1 cndtn1:g1 cndtn1:x1 grp1:1
## condition1 0.902
## group1 0.248 0.239
## exp1 0.426 0.453 0.143
## cndtn1:grp1 0.242 0.256 0.887 0.165
## condtn1:xp1 0.356 0.419 0.177 0.823 0.197
## group1:exp1 0.148 0.161 0.281 0.351 0.319 0.366
## cndtn1:g1:1 0.188 0.200 0.167 0.395 0.225 0.476 0.771
## optimizer (bobyqa) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# 'isSingular'
# all "REML = FALSE" to use ANOVA model comparison later on
# INITIAL model with indicator variables ---------------------------------
summary(m0_rateM <- glmer(rate ~ 1+
grpd.ungr+lhd.rhd+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_rd.rpt+lhd.rhd_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+grpd.ungr_rd.rpt|Subj)+
(1+lhd.rhd|item),
family="binomial",
control=glmerControl(optCtrl=list(maxfun=1e5),
optimizer="bobyqa"),
pDat_mm))## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: rate ~ 1 + grpd.ungr + lhd.rhd + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_rd.rpt + lhd.rhd_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## (1 + grpd.ungr + rd.rpt + grpd.ungr_rd.rpt | Subj) + (1 +
## lhd.rhd | item)
## Data: pDat_mm
## Control: glmerControl(optCtrl = list(maxfun = 1e+05), optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 1172.7 1281.1 -565.4 1130.7 1268
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5739 -0.1593 0.2867 0.4828 3.0638
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 2.57116 1.6035
## grpd.ungr 8.23428 2.8695 0.95
## rd.rpt 4.13984 2.0347 0.34 0.38
## grpd.ungr_rd.rpt 4.38362 2.0937 0.11 0.23 0.95
## item (Intercept) 0.01638 0.1280
## lhd.rhd 0.13368 0.3656 -1.00
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.5818 0.3629 4.359 1.31e-05 ***
## grpd.ungr 0.1798 0.6630 0.271 0.786291
## lhd.rhd 0.7818 0.6713 1.165 0.244157
## rd.rpt 2.3320 0.5704 4.089 4.34e-05 ***
## grpd.ungr_lhd.rhd 1.1029 1.2059 0.915 0.360408
## grpd.ungr_rd.rpt 3.5188 0.9340 3.767 0.000165 ***
## lhd.rhd_rd.rpt 1.4909 1.0070 1.481 0.138704
## grpd.ungr_lhd.rhd_rd.rpt 2.2744 1.5297 1.487 0.137071
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) grpd.n lhd.rh rd.rpt grpd.ngr_l. grpd.ngr_r. lhd._.
## grpd.ungr 0.902
## lhd.rhd 0.248 0.239
## rd.rpt 0.426 0.453 0.143
## grpd.ngr_l. 0.242 0.256 0.887 0.165
## grpd.ngr_r. 0.356 0.419 0.177 0.823 0.197
## lhd.rhd_rd. 0.148 0.161 0.281 0.351 0.319 0.366
## grpd.ng_._. 0.188 0.200 0.167 0.395 0.225 0.476 0.771
## optimizer (bobyqa) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# same model output
# two correlations close to 1 for Subj | -1 for items
summary(rePCA(m0_rateM))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4]
## Standard deviation 3.5300 2.5583 0.56875 0
## Proportion of Variance 0.6447 0.3386 0.01674 0
## Cumulative Proportion 0.6447 0.9833 1.00000 1
##
## $item
## Importance of components:
## [,1] [,2]
## Standard deviation 0.3874 3.472e-07
## Proportion of Variance 1.0000 0.000e+00
## Cumulative Proportion 1.0000 1.000e+00# 3 of 4 components for Subj | only intercept for item contribute (variance)
# try model reduction of randEff ---------------------------------
## start with ZERO CORRELATION model =================================
summary(zcp_rateM <- glmer(rate ~ 1+
grpd.ungr+lhd.rhd+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_rd.rpt+lhd.rhd_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+grpd.ungr_rd.rpt||Subj)+
(1+lhd.rhd||item),
family="binomial",
control=glmerControl(optCtrl=list(maxfun=1e5),
optimizer="bobyqa"),
pDat_mm))## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: rate ~ 1 + grpd.ungr + lhd.rhd + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_rd.rpt + lhd.rhd_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## (1 + grpd.ungr + rd.rpt + grpd.ungr_rd.rpt || Subj) + (1 +
## lhd.rhd || item)
## Data: pDat_mm
## Control: glmerControl(optCtrl = list(maxfun = 1e+05), optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 1202.9 1275.2 -587.5 1174.9 1275
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -4.2327 -0.2725 0.2721 0.4679 2.3564
##
## Random effects:
## Groups Name Variance Std.Dev.
## Subj (Intercept) 1.57842 1.2564
## Subj.1 grpd.ungr 4.03538 2.0088
## Subj.2 rd.rpt 2.28673 1.5122
## Subj.3 grpd.ungr_rd.rpt 1.75237 1.3238
## item (Intercept) 0.01059 0.1029
## item.1 lhd.rhd 0.18332 0.4282
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.4294 0.2607 5.484 4.16e-08 ***
## grpd.ungr -0.5452 0.4223 -1.291 0.197
## lhd.rhd 0.5521 0.5202 1.061 0.289
## rd.rpt 1.6939 0.3808 4.448 8.66e-06 ***
## grpd.ungr_lhd.rhd 0.7744 0.8570 0.904 0.366
## grpd.ungr_rd.rpt 2.1808 0.5363 4.066 4.78e-05 ***
## lhd.rhd_rd.rpt 1.1868 0.7656 1.550 0.121
## grpd.ungr_lhd.rhd_rd.rpt 0.8611 1.0715 0.804 0.422
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) grpd.n lhd.rh rd.rpt grpd.ngr_l. grpd.ngr_r. lhd._.
## grpd.ungr 0.045
## lhd.rhd 0.187 0.028
## rd.rpt 0.014 0.049 -0.020
## grpd.ngr_l. 0.036 0.184 0.041 0.049
## grpd.ngr_r. 0.098 0.032 0.051 0.175 0.031
## lhd.rhd_rd. 0.000 0.043 -0.033 0.289 0.047 0.118
## grpd.ng_._. 0.056 0.019 0.047 0.100 0.018 0.350 0.124 'Random effects:
Groups Name Variance Std.Dev.
Subj (Intercept) 1.57842 1.2564
Subj.1 grpd.ungr 4.03538 2.0088
Subj.2 rd.rpt 2.28673 1.5122
Subj.3 grpd.ungr_rd.rpt 1.75237 1.3238 # mean var. almost equal to sd
item (Intercept) 0.01059 0.1029 # almost no variance
item.1 lhd.rhd 0.18332 0.4282 # sd larger than variance
Number of obs: 1289, groups: Subj, 32; item, 24'## [1] "Random effects:\n Groups Name Variance Std.Dev.\n Subj (Intercept) 1.57842 1.2564 \n Subj.1 grpd.ungr 4.03538 2.0088 \n Subj.2 rd.rpt 2.28673 1.5122 \n Subj.3 grpd.ungr_rd.rpt 1.75237 1.3238 # mean var. almost equal to sd \n item (Intercept) 0.01059 0.1029 # almost no variance \n item.1 lhd.rhd 0.18332 0.4282 # sd larger than variance \n Number of obs: 1289, groups: Subj, 32; item, 24"## ANOVA model comparison - initial model - zcp model =================================
anova(m0_rateM, zcp_rateM)# m0_rateM has signif. better fit - removing all correlations too harsh
## back to initial model =================================
# try exclusion of randEff slopes with regards to zcp
### FIRST #############################
# try exclusion of "lhd.rhd" as slope per item
summary(m1_rateM <- glmer(rate ~ 1+
grpd.ungr+lhd.rhd+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_rd.rpt+lhd.rhd_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+grpd.ungr_rd.rpt|Subj)+
(1|item),
family="binomial",
control=glmerControl(optCtrl=list(maxfun=1e5),
optimizer="bobyqa"),
pDat_mm))## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: rate ~ 1 + grpd.ungr + lhd.rhd + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_rd.rpt + lhd.rhd_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## (1 + grpd.ungr + rd.rpt + grpd.ungr_rd.rpt | Subj) + (1 | item)
## Data: pDat_mm
## Control: glmerControl(optCtrl = list(maxfun = 1e+05), optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 1170.2 1268.3 -566.1 1132.2 1270
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.6493 -0.1576 0.2865 0.4938 3.1169
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 2.53998 1.5937
## grpd.ungr 8.12449 2.8503 0.95
## rd.rpt 4.06650 2.0166 0.34 0.38
## grpd.ungr_rd.rpt 4.31282 2.0767 0.12 0.24 0.95
## item (Intercept) 0.01979 0.1407
## Number of obs: 1289, groups: Subj, 32; item, 24
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.5740 0.3613 4.356 1.32e-05 ***
## grpd.ungr 0.1844 0.6599 0.279 0.779904
## lhd.rhd 0.7963 0.6634 1.200 0.229971
## rd.rpt 2.3249 0.5680 4.093 4.26e-05 ***
## grpd.ungr_lhd.rhd 1.0842 1.1897 0.911 0.362133
## grpd.ungr_rd.rpt 3.5059 0.9319 3.762 0.000168 ***
## lhd.rhd_rd.rpt 1.5060 0.9897 1.522 0.128083
## grpd.ungr_lhd.rhd_rd.rpt 2.2923 1.4935 1.535 0.124822
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) grpd.n lhd.rh rd.rpt grpd.ngr_l. grpd.ngr_r. lhd._.
## grpd.ungr 0.901
## lhd.rhd 0.260 0.241
## rd.rpt 0.430 0.456 0.147
## grpd.ngr_l. 0.244 0.270 0.900 0.169
## grpd.ngr_r. 0.359 0.421 0.181 0.822 0.200
## lhd.rhd_rd. 0.153 0.166 0.291 0.372 0.329 0.373
## grpd.ng_._. 0.195 0.206 0.176 0.405 0.235 0.510 0.796
## optimizer (bobyqa) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
summary(rePCA(m1_rateM))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4]
## Standard deviation 3.5108 2.5318 0.55486 0
## Proportion of Variance 0.6472 0.3366 0.01617 0
## Cumulative Proportion 0.6472 0.9838 1.00000 1
##
## $item
## Importance of components:
## [,1]
## Standard deviation 0.1407
## Proportion of Variance 1.0000
## Cumulative Proportion 1.0000# 1 of 4 components for Subj has zero variance | intercept only for item
### ANOVA model comparison - initial model - reduced model #############################
anova(m0_rateM, m1_rateM)# not significant, AIC smaller for m1_rateM - preferred
### SECOND #############################
# try exclusion of randEff for item altogether
summary(m2_rateM <- glmer(rate ~ 1+
grpd.ungr+lhd.rhd+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_rd.rpt+lhd.rhd_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt+grpd.ungr_rd.rpt|Subj),
family="binomial",
control=glmerControl(optCtrl=list(maxfun=1e5),
optimizer="bobyqa"),
pDat_mm))## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: rate ~ 1 + grpd.ungr + lhd.rhd + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_rd.rpt + lhd.rhd_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## (1 + grpd.ungr + rd.rpt + grpd.ungr_rd.rpt | Subj)
## Data: pDat_mm
## Control: glmerControl(optCtrl = list(maxfun = 1e+05), optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 1168.4 1261.3 -566.2 1132.4 1271
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.7248 -0.1564 0.2837 0.4990 3.1065
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 2.523 1.588
## grpd.ungr 8.062 2.839 0.95
## rd.rpt 4.045 2.011 0.34 0.38
## grpd.ungr_rd.rpt 4.295 2.073 0.12 0.24 0.95
## Number of obs: 1289, groups: Subj, 32
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.5686 0.3588 4.371 1.24e-05 ***
## grpd.ungr 0.1852 0.6553 0.283 0.777463
## lhd.rhd 0.7943 0.6614 1.201 0.229802
## rd.rpt 2.3150 0.5632 4.110 3.95e-05 ***
## grpd.ungr_lhd.rhd 1.0815 1.1858 0.912 0.361757
## grpd.ungr_rd.rpt 3.4936 0.9221 3.789 0.000151 ***
## lhd.rhd_rd.rpt 1.4963 0.9872 1.516 0.129608
## grpd.ungr_lhd.rhd_rd.rpt 2.2828 1.4906 1.531 0.125661
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) grpd.n lhd.rh rd.rpt grpd.ngr_l. grpd.ngr_r. lhd._.
## grpd.ungr 0.908
## lhd.rhd 0.261 0.242
## rd.rpt 0.433 0.461 0.148
## grpd.ngr_l. 0.245 0.271 0.901 0.170
## grpd.ngr_r. 0.363 0.427 0.183 0.832 0.202
## lhd.rhd_rd. 0.153 0.166 0.291 0.374 0.330 0.375
## grpd.ng_._. 0.195 0.207 0.177 0.407 0.236 0.514 0.796
## optimizer (bobyqa) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
summary(rePCA(m2_rateM))## $Subj
## Importance of components:
## [,1] [,2] [,3] [,4]
## Standard deviation 3.4996 2.5243 0.55304 0
## Proportion of Variance 0.6471 0.3367 0.01616 0
## Cumulative Proportion 0.6471 0.9838 1.00000 1# 1 of 4 components for Subj has zero variance
### ANOVA model comparison - m2_rateM - m1_rateM #############################
anova(m1_rateM, m2_rateM)# not significant, AIC smaller for m2_rateM - preferred
### THIRD #############################
# try exclusion of "rd.rpt" for Subj
summary(m3_rateM <- glmer(rate ~ 1+
grpd.ungr+lhd.rhd+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_rd.rpt+lhd.rhd_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+
(1+grpd.ungr+grpd.ungr_rd.rpt|Subj),
family="binomial",
control=glmerControl(optCtrl=list(maxfun=1e5),
optimizer="bobyqa"),
pDat_mm))## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: rate ~ 1 + grpd.ungr + lhd.rhd + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_rd.rpt + lhd.rhd_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## (1 + grpd.ungr + grpd.ungr_rd.rpt | Subj)
## Data: pDat_mm
## Control: glmerControl(optCtrl = list(maxfun = 1e+05), optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 1204.2 1276.5 -588.1 1176.2 1275
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -8.9531 -0.2444 0.3198 0.5289 2.9811
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 1.962 1.401
## grpd.ungr 5.799 2.408 0.91
## grpd.ungr_rd.rpt 3.459 1.860 -0.53 -0.22
## Number of obs: 1289, groups: Subj, 32
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.3083 0.2767 4.728 2.27e-06 ***
## grpd.ungr -0.3525 0.4953 -0.712 0.476644
## lhd.rhd 0.4416 0.5474 0.807 0.419790
## rd.rpt 1.8220 0.2307 7.899 2.82e-15 ***
## grpd.ungr_lhd.rhd 0.5074 0.9728 0.522 0.601954
## grpd.ungr_rd.rpt 2.2157 0.5788 3.828 0.000129 ***
## lhd.rhd_rd.rpt 1.3443 0.4454 3.018 0.002543 **
## grpd.ungr_lhd.rhd_rd.rpt 1.9777 1.1119 1.779 0.075286 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) grpd.n lhd.rh rd.rpt grpd.ngr_l. grpd.ngr_r. lhd._.
## grpd.ungr 0.822
## lhd.rhd 0.167 0.148
## rd.rpt 0.069 0.025 0.019
## grpd.ngr_l. 0.145 0.193 0.835 0.026
## grpd.ngr_r. -0.239 -0.033 0.007 0.370 0.022
## lhd.rhd_rd. 0.045 0.046 0.028 0.542 0.012 0.336
## grpd.ng_._. 0.004 0.022 -0.279 0.317 -0.097 0.394 0.399summary(rePCA(m3_rateM))## $Subj
## Importance of components:
## [,1] [,2] [,3]
## Standard deviation 2.8351 1.7611 0.28357
## Proportion of Variance 0.7164 0.2764 0.00717
## Cumulative Proportion 0.7164 0.9928 1.00000# all components seem okay
### ANOVA model comparison - m3_rateM - m2_rateM #############################
anova(m2_rateM, m3_rateM)# significant, AIC smaller for m2_rateM - preferred
### FOURTH #############################
# back to m2_rateM
# try exclusion of "grpd.ungr_rd.rpt" for Subj
summary(m4_rateM <- glmer(rate ~ 1+
grpd.ungr+lhd.rhd+rd.rpt+
grpd.ungr_lhd.rhd+grpd.ungr_rd.rpt+lhd.rhd_rd.rpt+
grpd.ungr_lhd.rhd_rd.rpt+
(1+grpd.ungr+rd.rpt|Subj),
family="binomial",
control=glmerControl(optCtrl=list(maxfun=1e5),
optimizer="bobyqa"),
pDat_mm))## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: rate ~ 1 + grpd.ungr + lhd.rhd + rd.rpt + grpd.ungr_lhd.rhd +
## grpd.ungr_rd.rpt + lhd.rhd_rd.rpt + grpd.ungr_lhd.rhd_rd.rpt +
## (1 + grpd.ungr + rd.rpt | Subj)
## Data: pDat_mm
## Control: glmerControl(optCtrl = list(maxfun = 1e+05), optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 1175.8 1248.1 -573.9 1147.8 1275
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -4.9809 -0.1912 0.2966 0.4746 3.0101
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 2.098 1.449
## grpd.ungr 6.509 2.551 0.93
## rd.rpt 2.417 1.555 0.27 0.34
## Number of obs: 1289, groups: Subj, 32
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.4875 0.2950 5.042 4.60e-07 ***
## grpd.ungr -0.1468 0.5171 -0.284 0.777
## lhd.rhd 0.6953 0.5797 1.199 0.230
## rd.rpt 1.9285 0.3983 4.842 1.29e-06 ***
## grpd.ungr_lhd.rhd 0.8638 1.0209 0.846 0.397
## grpd.ungr_rd.rpt 2.9187 0.5246 5.563 2.65e-08 ***
## lhd.rhd_rd.rpt 1.3290 0.7716 1.722 0.085 .
## grpd.ungr_lhd.rhd_rd.rpt 2.0324 0.9946 2.043 0.041 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) grpd.n lhd.rh rd.rpt grpd.ngr_l. grpd.ngr_r. lhd._.
## grpd.ungr 0.835
## lhd.rhd 0.188 0.146
## rd.rpt 0.247 0.291 0.051
## grpd.ngr_l. 0.147 0.181 0.842 0.083
## grpd.ngr_r. 0.138 0.108 0.078 0.382 0.071
## lhd.rhd_rd. 0.065 0.088 0.194 0.350 0.274 0.306
## grpd.ng_._. 0.107 0.101 0.092 0.326 0.074 0.600 0.374summary(rePCA(m4_rateM))## $Subj
## Importance of components:
## [,1] [,2] [,3]
## Standard deviation 2.9565 1.4417 0.4529
## Proportion of Variance 0.7929 0.1885 0.0186
## Cumulative Proportion 0.7929 0.9814 1.0000# all components seem okay
### ANOVA model comparison - m4_rateM - m2_rateM #############################
anova(m4_rateM, m2_rateM)# significant, AIC smaller for m2_rateM - preferred
# no better model fit achievable
# use m2_rateM
# FINAL model ---------------------------------
## model with factors and * for fixed effects =================================
# use REML=TRUE
summary(rateM <- glmer(rate ~ condition*group*exp+
(1+grpd.ungr+rd.rpt+grpd.ungr_rd.rpt|Subj),
family="binomial",
control=glmerControl(optCtrl=list(maxfun=1e5),
optimizer="bobyqa"),
pDat_mm))## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula:
## rate ~ condition * group * exp + (1 + grpd.ungr + rd.rpt + grpd.ungr_rd.rpt |
## Subj)
## Data: pDat_mm
## Control: glmerControl(optCtrl = list(maxfun = 1e+05), optimizer = "bobyqa")
##
## AIC BIC logLik deviance df.resid
## 1168.4 1261.3 -566.2 1132.4 1271
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.7248 -0.1564 0.2837 0.4990 3.1065
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Subj (Intercept) 2.523 1.588
## grpd.ungr 8.062 2.839 0.95
## rd.rpt 4.045 2.011 0.34 0.38
## grpd.ungr_rd.rpt 4.295 2.073 0.12 0.24 0.95
## Number of obs: 1289, groups: Subj, 32
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.5686 0.3588 4.371 1.24e-05 ***
## condition1 0.1852 0.6553 0.283 0.777463
## group1 0.7943 0.6614 1.201 0.229802
## exp1 2.3150 0.5632 4.110 3.95e-05 ***
## condition1:group1 1.0815 1.1858 0.912 0.361757
## condition1:exp1 3.4936 0.9221 3.789 0.000151 ***
## group1:exp1 1.4963 0.9872 1.516 0.129608
## condition1:group1:exp1 2.2828 1.4906 1.531 0.125661
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) cndtn1 group1 exp1 cndtn1:g1 cndtn1:x1 grp1:1
## condition1 0.908
## group1 0.261 0.242
## exp1 0.433 0.461 0.148
## cndtn1:grp1 0.245 0.271 0.901 0.170
## condtn1:xp1 0.363 0.427 0.183 0.832 0.202
## group1:exp1 0.153 0.166 0.291 0.374 0.330 0.375
## cndtn1:g1:1 0.195 0.207 0.177 0.407 0.236 0.514 0.796
## optimizer (bobyqa) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')# isSingular
save(rateM, file = "rateM.Rdata")
# *Nice table output* ---------------------------------
# sjPlot::tab_model((rateM), show.se = TRUE, show.stat = TRUE, collapse.ci=TRUE, transform = NULL, file="rateM.doc")
# clean workspace ---------------------------------
# ls()
rm(m0_rateM, m1_rateM, m2_rateM, m3_rateM, m4_rateM, zcp_rateM)9 Production accuracy model
9.1 Visualizations and Post-Hoc-Tests
# interaction visualization ---------------------------------
load("rateM.Rdata")
# condition:exp =================================
### condition ~ exp #############################
emmeans::emmip(rateM, condition ~ exp, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - production accuracy",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
### exp ~ condition #############################
emmeans::emmip(rateM, exp ~ condition, CIs = TRUE,
CIarg = list(lwd = 1.1, alpha=0.7),
linearg = list(linetype = "dashed", alpha=0)) +
labs(
y = "Model predictions - production accuracy",
) +
scale_color_brewer(palette = "Set1") +
theme(axis.text=element_text(size=11),
axis.title.x=element_text(size=11),
axis.title.y=element_text(size=11),
legend.justification = c("right", "top"),
legend.title = element_text(size=11, face="bold"),
legend.text=element_text(size=11, face="bold")
)
# POST-HOC testing ---------------------------------
# what are the significant differences
emmeans::emmeans(rateM, pairwise ~ condition*exp) %>%
summary(by = NULL, adjust = "bonferroni")## $emmeans
## condition exp emmean SE df asymp.LCL asymp.UCL
## grouped reading_aloud 3.69 0.994 Inf 1.210 6.17
## ungrouped reading_aloud 1.76 0.232 Inf 1.181 2.34
## grouped repetition -0.37 0.627 Inf -1.936 1.20
## ungrouped repetition 1.19 0.200 Inf 0.693 1.69
##
## Results are averaged over the levels of: group
## Results are given on the logit (not the response) scale.
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 4 estimates
##
## $contrasts
## contrast estimate SE df z.ratio
## grouped reading_aloud - ungrouped reading_aloud 1.932 0.949 Inf 2.037
## grouped reading_aloud - grouped repetition 4.062 0.981 Inf 4.141
## grouped reading_aloud - ungrouped repetition 2.500 1.043 Inf 2.398
## ungrouped reading_aloud - grouped repetition 2.130 0.637 Inf 3.342
## ungrouped reading_aloud - ungrouped repetition 0.568 0.312 Inf 1.820
## grouped repetition - ungrouped repetition -1.562 0.620 Inf -2.520
## p.value
## 0.2501
## 0.0002
## 0.0989
## 0.0050
## 0.4124
## 0.0705
##
## Results are averaged over the levels of: group
## Results are given on the log odds ratio (not the response) scale.
## P value adjustment: bonferroni method for 6 tests10 Raw data plots
10.1 Production accuracy
10.2 Plotting means and error bars
Functions for mean + CI plots
+ start with this script chunk “SEsummary” + run it and then go on to
the next chunk to plot the data on all cues “cues separately”
# summarySE between ---------------------------------
## Gives count, mean, standard deviation, standard error of the mean, and confidence interval (default 95%)
## data: a data frame
## measurevar: the name of a column that contains the variable to be summarized
## groupvars: a vector containing names of columns that contain grouping variables
## na.rm: a boolean that indicates whether to ignore NA's
## conf.interval: the percent range of the confidence interval (default is 95%)
summarySE <- function(data=NULL, measurevar, groupvars=NULL, na.rm=TRUE,
conf.interval=.95, .drop=TRUE) {
# New version of length which can handle NA's: if na.rm==T, don't count them
length2 <- function (x, na.rm=TRUE) {
if (na.rm) sum(!is.na(x))
else length(x)
}
# This does the summary. For each group's data frame, return a vector with
# N, mean, and sd
datac <- ddply(data, groupvars, .drop=.drop,
.fun = function(xx, col) {
c(N = length2(xx[[col]], na.rm=na.rm),
mean = mean (xx[[col]], na.rm=na.rm),
sd = sd (xx[[col]], na.rm=na.rm)
)
},
measurevar
)
# Rename the "mean" column
datac <- plyr::rename(datac, c("mean" = measurevar))
datac$se <- datac$sd / sqrt(datac$N) # Calculate standard error of the mean
# Confidence interval multiplier for standard error
# Calculate t-statistic for confidence interval:
# e.g., if conf.interval is .95, use .975 (above/below), and use df=N-1
ciMult <- qt(conf.interval/2 + .5, datac$N-1)
datac$ci <- datac$se * ciMult
return(datac)
}
# normDataWithin ---------------------------------
## Norms the data within specified groups in a data frame; it normalizes each
## subject (identified by idvar) so that they have the same mean, within each group
## specified by betweenvars
## data: a data frame.
## idvar: the name of a column that identifies each subject (or matched subjects)
## measurevar: the name of a column that contains the variable to be summariezed
## betweenvars: a vector containing names of columns that are between-subjects variables
## na.rm: a boolean that indicates whether to ignore NA's
normDataWithin <- function(data=NULL, idvar, measurevar, betweenvars=NULL,
na.rm=TRUE, .drop=TRUE) {
# Measure var on left, idvar + between vars on right of formula
data.subjMean <- ddply(data, c(idvar, betweenvars), .drop=.drop,
.fun = function(xx, col, na.rm) {
c(subjMean = mean(xx[,col], na.rm=na.rm))
},
measurevar,
na.rm
)
# Put the subject means with original data
data <- merge(data, data.subjMean)
# Get the normalized data in a new column
measureNormedVar <- paste(measurevar, "_norm", sep="")
data[,measureNormedVar] <- data[,measurevar] - data[,"subjMean"] +
mean(data[,measurevar], na.rm=na.rm)
# Remove this subject mean column
data$subjMean <- NULL
return(data)
}
# SEsummarywithin ---------------------------------
## Summarizes data, handling within-subjects variables by removing inter-subject variability.
## It will still work if there are no within-S variables.
## Gives count, un-normed mean, normed mean (with same between-group mean),
## standard deviation, standard error of the mean, and confidence interval.
## If there are within-subject variables, calculate adjusted values using method from Morey (2008)
## data: a data frame.
## measurevar: the name of a column that contains the variable to be summariezed
## betweenvars: a vector containing names of columns that are between-subjects variables
## withinvars: a vector containing names of columns that are within-subjects variables
## idvar: the name of a column that identifies each subject (or matched subjects)
## na.rm: a boolean that indicates whether to ignore NA's
## conf.interval: the percent range of the confidence interval (default is 95%)
summarySEwithin <- function(data=NULL, measurevar, betweenvars=NULL, withinvars=NULL,
idvar=NULL, na.rm=TRUE, conf.interval=.95, .drop=TRUE) {
# Ensure that the betweenvars and withinvars are factors
factorvars <- vapply(data[, c(betweenvars, withinvars), drop=FALSE],
FUN=is.factor, FUN.VALUE=logical(1))
if (!all(factorvars)) {
nonfactorvars <- names(factorvars)[!factorvars]
message("Automatically converting the following non-factors to factors: ",
paste(nonfactorvars, collapse = ", "))
data[nonfactorvars] <- lapply(data[nonfactorvars], factor)
}
# Get the means from the un-normed data
datac <- summarySE(data, measurevar, groupvars=c(betweenvars, withinvars),
na.rm=na.rm, conf.interval=conf.interval, .drop=.drop)
# Drop all the unused columns (these will be calculated with normed data)
datac$sd <- NULL
datac$se <- NULL
datac$ci <- NULL
# Norm each subject's data
ndata <- normDataWithin(data, idvar, measurevar, betweenvars, na.rm, .drop=.drop)
# This is the name of the new column
measurevar_n <- paste(measurevar, "_norm", sep="")
# Collapse the normed data - now we can treat between and within vars the same
ndatac <- summarySE(ndata, measurevar_n, groupvars=c(betweenvars, withinvars),
na.rm=na.rm, conf.interval=conf.interval, .drop=.drop)
# Apply correction from Morey (2008) to the standard error and confidence interval
# Get the product of the number of conditions of within-S variables
nWithinGroups <- prod(vapply(ndatac[,withinvars, drop=FALSE], FUN=nlevels,
FUN.VALUE=numeric(1)))
correctionFactor <- sqrt( nWithinGroups / (nWithinGroups-1) )
# Apply the correction factor
ndatac$sd <- ndatac$sd * correctionFactor
ndatac$se <- ndatac$se * correctionFactor
ndatac$ci <- ndatac$ci * correctionFactor
# Combine the un-normed means with the normed results
merge(datac, ndatac)
}Mean an CI plot for cues on name2
# Cues ---------------------------------
Cues <- pDat %>%
dplyr::select(Subj, item, group, exp, condition, accuracy, f0_range2, s8reln2, pause3) %>%
tidyr::gather("variable", "measurement", -Subj, -item, -group, -exp, -condition, -accuracy)
Cues$group <- as.factor(as.character(Cues$group))
Cues$exp <- as.factor(as.character(Cues$exp))
Cues$condition <- as.factor(as.character(Cues$condition))
Cues$accuracy <- as.factor(as.character(Cues$accuracy))
Cues$variable <- as.factor(as.character(Cues$variable))
Cues$measurement <- as.numeric(as.character(Cues$measurement))
## f0_range =================================
f0.range <- Cues %>%
filter(variable=="f0_range2") %>%
group_by(Subj, group, exp, condition, accuracy, variable) %>%
summarySEwithin(measurevar = "measurement",
betweenvars = "group",
withinvars = c("variable", "condition", "exp", "accuracy"),
idvar = "Subj") %>%
ggplot(aes(x=accuracy,y=measurement, shape = condition, color=condition)) +
geom_point(position = position_dodge2(width = .5)) +
facet_grid(exp~group)+
geom_errorbar(aes(ymax = measurement + ci, ymin = measurement - ci), position = position_dodge2(width = .2), width=.5) +
theme_bw() +
theme(axis.text=element_text(size=14),
axis.title=element_text(size=14),
strip.text.x = element_text(size=14),
strip.text.y = element_blank(),
legend.title=element_text(size=14),
legend.text=element_text(size=14)) +
labs(
x = "",
y = "f0-range on name2 (st)"
)
f0.range
# ggsave("f0_range2.jpg", width = 18, height = 14, dpi = 300, units = "cm")
## final lengthening =================================
fin.length <- Cues %>%
filter(variable=="s8reln2") %>%
group_by(Subj, group, exp, condition, accuracy, variable) %>%
summarySEwithin(measurevar = "measurement",
betweenvars = "group",
withinvars = c("variable", "condition", "exp", "accuracy"),
idvar = "Subj") %>%
ggplot(aes(x=accuracy,y=measurement, shape = condition, color=condition)) +
geom_point(position = position_dodge2(width = .5)) +
facet_grid(exp~group)+
geom_errorbar(aes(ymax = measurement + ci, ymin = measurement - ci), position = position_dodge2(width = .2), width=.5) +
theme_bw() +
theme(axis.text=element_text(size=14),
axis.title=element_text(size=14),
strip.text.x = element_text(size=14),
strip.text.y = element_blank(),
legend.title=element_text(size=14),
legend.text=element_text(size=14)) +
labs(
x = "",
y = "final lengthening on name2 (%)"
)
fin.length
# ggsave("fin-length.jpg", width = 18, height = 14, dpi = 300, units = "cm")
## pause 3 =================================
pause3 <- Cues %>%
filter(variable=="pause3") %>%
group_by(Subj, group, exp, condition, accuracy, variable) %>%
summarySEwithin(measurevar = "measurement",
betweenvars = "group",
withinvars = c("variable", "condition", "exp", "accuracy"),
idvar = "Subj") %>%
ggplot(aes(x=accuracy,y=measurement, shape = condition, color=condition)) +
geom_point(position = position_dodge2(width = .5)) +
facet_grid(exp~group)+
geom_errorbar(aes(ymax = measurement + ci, ymin = measurement - ci), position = position_dodge2(width = .2), width=.5) +
theme_bw() +
theme(axis.text=element_text(size=14),
axis.title=element_text(size=14),
strip.text.x = element_text(size=14),
strip.text.y = element_text(size=14),
legend.title=element_text(size=14),
legend.text=element_text(size=14)) +
labs(
x = "",
y = "pause after name2 (%)"
)
pause3
# ggsave("pause3.jpg", width = 18, height = 14, dpi = 300, units = "cm")
## plot grids =================================
library(ggpubr)
pDat_cues <- ggarrange(f0.range,fin.length,pause3, ncol=3, nrow=1, common.legend = TRUE, legend="bottom")
pDat_cues
# ggsave("cue_grid.jpg", plot = pDat_cues, width = 18, height = 10, dpi = 300)
## clean workspace =================================
rm(Cues, f0.range, fin.length, pause3, pDat_cues, normDataWithin, summarySE, summarySEwithin)11 Correlations - working memory measures / accuracy
# import working memory data ---------------------------------
workMem <- read.csv2("data/workMem.csv", header = TRUE)
# enumerate accuracy variable ---------------------------------
pDat$rating_match <- as.numeric(as.factor(pDat$rating_match)) -1
## participant check =================================
## list first entry for exp for each participant
pDat %>%
group_by(Subj, group, exp) %>%
dplyr::summarise(mean = round(mean(rating_match)/11*100, digits=2)) %>%
group_by(Subj) %>%
summarise(exp = first(exp)) %>%
as.data.frame()### RHDP11 only took part in repetition exp
## list how many entries for exp for each participant exist
pDat %>%
group_by(Subj, group, exp) %>%
dplyr::summarise(mean = round(mean(rating_match)/11*100, digits=2)) %>%
group_by(Subj) %>%
summarise(exp = n_distinct(exp)) %>%
as.data.frame()# LHDP01, LHDP02, LHDP08, LHDP12, RHDP19, RHDP20 only took part in reading_aloud
## working memory subset =================================
## only individuals that participated in both exp
workMemBoth <- subset(workMem,
workMem$Subj=="LHDP03" | workMem$Subj=="LHDP04" | workMem$Subj=="LHDP07" | workMem$Subj=="LHDP09" |
workMem$Subj=="LHDP10" | workMem$Subj=="LHDP11" | workMem$Subj=="LHDP13" | workMem$Subj=="LHDP14" |
workMem$Subj=="LHDP15" | workMem$Subj=="LHDP16" | workMem$Subj=="RHDP01" | workMem$Subj=="RHDP02" |
workMem$Subj=="RHDP03" | workMem$Subj=="RHDP04" | workMem$Subj=="RHDP05" | workMem$Subj=="RHDP08" |
workMem$Subj=="RHDP09" | workMem$Subj=="RHDP12" | workMem$Subj=="RHDP13" | workMem$Subj=="RHDP14" |
workMem$Subj=="RHDP15" | workMem$Subj=="RHDP16" | workMem$Subj=="RHDP17" | workMem$Subj=="RHDP18" |
workMem$Subj=="RHDP22")
## data subset =================================
## only individuals that participated in both exp
meanBoth <- pDat %>%
group_by(Subj) %>%
filter(n_distinct(exp)==2)
# data aggregation ---------------------------------
### means per Subj, group, exp, condition merged with working memory measures #############################
mean_rating_Subj <- meanBoth %>%
group_by(Subj, group, exp, condition) %>%
dplyr::summarise(mean = round(mean(rating_match)/11*100, digits=2))
allMeans <- merge(mean_rating_Subj, workMemBoth[, c("Subj", "compDigit")], by=c("Subj"), all.x=TRUE)
### means per Subj, group, exp - overall conditions merged with working memory measures #############################
mean_rating_Subj1 <- meanBoth %>%
group_by(Subj, group, exp) %>%
dplyr::summarise(mean = round(mean(rating_match)/11*100, digits=2))
allMeans1 <- merge(mean_rating_Subj1, workMemBoth[, c("Subj", "compDigit")], by=c("Subj"), all.x=TRUE)
# plots ---------------------------------
### compDigit - COMPOSITE SCORE DIGIT SPAN forward & backward #############################
## REPETITION =================================
# mean accuracy assessment over group and OVERALL conditions
### these correlational measures are REPORTED IN PAPER ("3.2.3 Effects of Working Memory") ###
corrRep <- facet(ggscatter(subset(allMeans1, allMeans1$exp=="repetition"), x="mean", y="compDigit",
add = "reg.line", conf.int = TRUE,
cor.coef = TRUE, cor.coeff.args = list(size=6, label.x = 50, label.y = 80, label.sep = "\n"), cor.method = "spearman",
xlab = "mean accuracy assessment repetition overall", ylab = "composite score digit span")+
theme_bw(),
facet.by = c("group"),
short.panel.labs = TRUE, # Allow long labels in panels
panel.labs.background = list(fill = "lightgrey", color = "lightgrey"),
panel.labs.font = list(face = NULL, color = "black", size = 20, angle = NULL))+
font("xlab", size = 20, color = "black")+
font("ylab", size = 20, color = "black")+
font("xy.text", size = 19, color = "black")
corrRep
# ggsave("correlations_repetition_overall.jpg", plot = corrRep1, width = 18, height = 16, dpi = 300)
# mean accuracy assessment over condition AND group
corrRep1 <- facet(ggscatter(subset(allMeans, allMeans$exp=="repetition"), x="mean", y="compDigit",
add = "reg.line", conf.int = TRUE,
cor.coef = TRUE, cor.coeff.args = list(size=6, label.x = 25, label.y = 75, label.sep = "\n"), cor.method = "spearman",
xlab = "mean accuracy assessment repetition", ylab = "composite score digit span")+
theme_bw(),
facet.by = c("group","condition"),
short.panel.labs = TRUE, # Allow long labels in panels
panel.labs.background = list(fill = "lightgrey", color = "lightgrey"),
panel.labs.font = list(face = NULL, color = "black", size = 20, angle = NULL))+
font("xlab", size = 20, color = "black")+
font("ylab", size = 20, color = "black")+
font("xy.text", size = 19, color = "black")
corrRep1
# ggsave("correlations_repetition.jpg", plot = corrRep, width = 18, height = 16, dpi = 300)
# READING ALOUD =================================
# mean accuracy assessment over group and OVERALL conditions
corrRead <- facet(ggscatter(subset(allMeans1, allMeans1$exp=="reading_aloud"), x="mean", y="compDigit",
add = "reg.line", conf.int = TRUE,
cor.coef = TRUE, cor.coeff.args = list(size=6, label.x = 50, label.y = 80, label.sep = "\n"), cor.method = "spearman",
xlab = "mean accuracy assessment reading_aloud overall", ylab = "composite score digit span")+
theme_bw(),
facet.by = c("group"),
short.panel.labs = TRUE, # Allow long labels in panels
panel.labs.background = list(fill = "lightgrey", color = "lightgrey"),
panel.labs.font = list(face = NULL, color = "black", size = 20, angle = NULL))+
font("xlab", size = 20, color = "black")+
font("ylab", size = 20, color = "black")+
font("xy.text", size = 19, color = "black")
corrRead
# ggsave("correlations_reading_aloud_overall.jpg", plot = corrRead1, width = 18, height = 16, dpi = 300)
# mean accuracy assessment over condition AND group
corrRead1 <- facet(ggscatter(subset(allMeans, allMeans$exp=="reading_aloud"), x="mean", y="compDigit",
add = "reg.line", conf.int = TRUE,
cor.coef = TRUE, cor.coeff.args = list(size=6, label.x = 20, label.y = 75, label.sep = "\n"), cor.method = "spearman",
xlab = "mean accuracy assessment reading_aloud", ylab = "composite score digit span")+
theme_bw(),
facet.by = c("group","condition"),
short.panel.labs = TRUE, # Allow long labels in panels
panel.labs.background = list(fill = "lightgrey", color = "lightgrey"),
panel.labs.font = list(face = NULL, color = "black", size = 20, angle = NULL))+
font("xlab", size = 20, color = "black")+
font("ylab", size = 20, color = "black")+
font("xy.text", size = 19, color = "black")
corrRead1
# ggsave("correlations_reading_aloud.jpg", plot = corrRead, width = 18, height = 16, dpi = 300)12 Summary Tables of Cue Use for Patient Groups and Model-Stimuli
require(flextable)
data_corr <- subset(data, data$accuracy=="correct")
# generate descriptive summary for pause 3 per group (pause3 duration relative to utterance duration in %) ---------------------------------
pause3group <- data_corr %>%
dplyr::group_by(exp, group, condition) %>%
dplyr::summarise(mean = round(mean(pause3),3), n = n(), min = round(min(pause3),3),
max = round(max(pause3),3), sd = round(sd(pause3),3))
## generate and safe respective table for pause 3 per group =================================
pause3groupFT <- flextable(pause3group) %>%
merge_v(j = c("exp", "group")) %>%
valign(valign = "top") %>%
theme_box() %>% set_table_properties(layout = "autofit")
pause3groupFT <- add_header_row(
x = pause3groupFT, values = c("Descriptive summary for pause3 by patient group and model stimuli
Only correct production accuracy"),
colwidths = c(8))
pause3groupFT <- align(pause3groupFT, i = 1, part = "header", align = "center")
pause3groupFT <- bold(pause3groupFT, j = 1, i = ~ !is.na(exp))
pause3groupFT <- fontsize(pause3groupFT, part = "all", size = 12)
pause3groupFT <- theme_vanilla(pause3groupFT)
pause3groupFTDescriptive summary for pause3 by patient group and model stimuli | |||||||
exp | group | condition | mean | n | min | max | sd |
reading_aloud | LHDP | grouped | 24.719 | 137 | 0.000 | 47.411 | 9.037 |
ungrouped | 2.718 | 130 | 0.000 | 18.251 | 4.330 | ||
RHDP | grouped | 23.150 | 145 | 0.000 | 44.024 | 8.801 | |
ungrouped | 0.855 | 160 | 0.000 | 11.233 | 2.445 | ||
repetition | LHDP | grouped | 16.952 | 56 | 5.429 | 36.583 | 7.759 |
ungrouped | 1.218 | 74 | 0.000 | 16.560 | 2.855 | ||
model_stimuli | grouped | 14.801 | 6 | 8.818 | 23.269 | 5.716 | |
ungrouped | 2.265 | 6 | 0.000 | 4.773 | 2.486 | ||
RHDP | grouped | 14.302 | 84 | 0.000 | 32.290 | 8.124 | |
ungrouped | 0.508 | 139 | 0.000 | 11.728 | 1.828 | ||
save_as_docx(
"Table Summary Pause Cue" = pause3groupFT,
path = "group_sum_pause.docx")
rm(pause3group)
rm(pause3groupFT)
# F0-range ---------------------------------
# generate descriptive summary for f0_range2 per group (difference in semitones btw L and H annotated on name2) ---------------------------------
f0_range2group <- data_corr %>%
dplyr::group_by(exp, group, condition) %>%
dplyr::summarise(mean = round(mean(f0_range2, na.rm=T),3), n = n(), miss = sum(is.na(f0_range2)), min = round(min(f0_range2, na.rm=T),3),
max = round(max(f0_range2, na.rm=T),3), sd = round(sd(f0_range2, na.rm=T),3))
## generate and safe respective table for f0_range2 per group =================================
f0_range2groupFT <- flextable(f0_range2group) %>%
merge_v(j = c("exp", "group")) %>%
valign(valign = "top") %>%
theme_box() %>% set_table_properties(layout = "autofit")
f0_range2groupFT <- add_header_row(
x = f0_range2groupFT, values = c("Descriptive summary for f0-range on name2 by patient group and model stimuli
Only correct production accuracy"),
colwidths = c(9))
f0_range2groupFT <- align(f0_range2groupFT, i = 1, part = "header", align = "center")
f0_range2groupFT <- bold(f0_range2groupFT, j = 1, i = ~ !is.na(exp))
f0_range2groupFT <- fontsize(f0_range2groupFT, part = "all", size = 12)
f0_range2groupFT <- theme_vanilla(f0_range2groupFT)
f0_range2groupFTDescriptive summary for f0-range on name2 by patient group and model stimuli | ||||||||
exp | group | condition | mean | n | miss | min | max | sd |
reading_aloud | LHDP | grouped | 10.975 | 137 | 4 | 3.830 | 18.96 | 3.103 |
ungrouped | 5.836 | 130 | 3 | -2.570 | 13.62 | 2.650 | ||
RHDP | grouped | 6.807 | 145 | 7 | -0.001 | 18.92 | 4.122 | |
ungrouped | 3.429 | 160 | 8 | -3.300 | 20.11 | 2.744 | ||
repetition | LHDP | grouped | 7.051 | 56 | 1 | 0.600 | 12.94 | 2.565 |
ungrouped | 2.464 | 74 | 0 | -3.280 | 7.23 | 1.722 | ||
model_stimuli | grouped | 9.825 | 6 | 0 | 8.920 | 10.96 | 0.676 | |
ungrouped | 4.953 | 6 | 0 | 4.240 | 5.76 | 0.630 | ||
RHDP | grouped | 6.990 | 84 | 2 | 1.140 | 19.24 | 3.734 | |
ungrouped | 2.930 | 139 | 10 | -10.300 | 15.25 | 3.208 | ||
save_as_docx(
"Table Summary f0-Range Cue" = f0_range2groupFT,
path = "group_sum_f0.docx")
rm(f0_range2group)
rm(f0_range2groupFT)
# Final lengthening ---------------------------------
# generate descriptive summary for s8reln2 (length of final vowel, s8, in % relative to duration of name2) ---------------------------------
s8reln2group <- data_corr %>%
dplyr::group_by(exp, group, condition) %>%
dplyr::summarise(mean = round(mean(s8reln2),3), n = n(), min = round(min(s8reln2),3),
max = round(max(s8reln2),3), sd = round(sd(s8reln2),3))
## 8B generate and safe respective table for s8reln2 =================================
s8reln2groupFT <- flextable(s8reln2group) %>%
merge_v(j = c("exp", "group")) %>%
valign(valign = "top") %>%
theme_box() %>% set_table_properties(layout = "autofit")
s8reln2groupFT <- add_header_row(
x = s8reln2groupFT, values = c("Descriptive summary for final lengthening on name2 by patient group and model stimuli
Only correct production accuracy"),
colwidths = c(8))
s8reln2groupFT <- align(s8reln2groupFT, i = 1, part = "header", align = "center")
s8reln2groupFT <- bold(s8reln2groupFT, j = 1, i = ~ !is.na(exp))
s8reln2groupFT <- fontsize(s8reln2groupFT, part = "all", size = 12)
s8reln2groupFT <- theme_vanilla(s8reln2groupFT)
s8reln2groupFTDescriptive summary for final lengthening on name2 by patient group and model stimuli | |||||||
exp | group | condition | mean | n | min | max | sd |
reading_aloud | LHDP | grouped | 49.336 | 137 | 32.896 | 65.693 | 7.067 |
ungrouped | 40.581 | 130 | 19.269 | 62.539 | 10.529 | ||
RHDP | grouped | 47.747 | 145 | 34.456 | 64.619 | 6.043 | |
ungrouped | 37.626 | 160 | 15.122 | 64.565 | 8.957 | ||
repetition | LHDP | grouped | 52.923 | 56 | 39.025 | 66.727 | 5.920 |
ungrouped | 38.019 | 74 | 17.605 | 59.238 | 9.600 | ||
model_stimuli | grouped | 42.328 | 6 | 36.938 | 48.057 | 4.334 | |
ungrouped | 34.232 | 6 | 26.134 | 40.834 | 5.958 | ||
RHDP | grouped | 46.885 | 84 | 33.667 | 65.794 | 6.054 | |
ungrouped | 35.003 | 139 | 20.271 | 53.724 | 7.714 | ||
save_as_docx(
"Table Summary Final Lengthening Cue" = s8reln2groupFT,
path = "group_sum_FL.docx")
rm(s8reln2group)
rm(s8reln2groupFT)13 Boxplot model-stimuli
mod_stim <- subset(data_corr, data_corr$Subj=="model_stimulus")
# cue value distribution by condition significantly different from each other?
mod_grou <- subset(mod_stim, mod_stim$condition=="grouped")
mod_ungrou <- subset(mod_stim, mod_stim$condition=="ungrouped")
# vectors / groups for cues
x1 <- mod_grou$f0_range2
x2 <- mod_ungrou$f0_range2
y1 <- mod_grou$s8reln2
y2 <- mod_ungrou$s8reln2
z1 <- mod_grou$pause3
z2 <- mod_ungrou$pause3
# dependent 2-group Wilcoxon Signed Rank Test
# where both groups are numeric
wilcox.test(x1,x2,paired=TRUE) ##
## Wilcoxon signed rank exact test
##
## data: x1 and x2
## V = 21, p-value = 0.03125
## alternative hypothesis: true location shift is not equal to 0wilcox.test(y1,y2,paired=TRUE)##
## Wilcoxon signed rank exact test
##
## data: y1 and y2
## V = 21, p-value = 0.03125
## alternative hypothesis: true location shift is not equal to 0wilcox.test(z1,z2,paired=TRUE)##
## Wilcoxon signed rank exact test
##
## data: z1 and z2
## V = 21, p-value = 0.03125
## alternative hypothesis: true location shift is not equal to 0# boxplot f0-range2
f0 <- ggplot(mod_stim, aes(condition, f0_range2))
f0_p <- f0 + geom_boxplot(outlier.shape = NA) +
geom_jitter(width = 0.1) +
labs(
title = "F0 range on name2 by condition for model stimuli",
) +
theme_bw() +
theme(axis.text.x = element_text(size=14),
axis.text.y = element_text(size=14),
axis.title=element_text(size=14)) +
annotate("segment", x = 1, xend = 2, y = 11.5, yend = 11.5,
colour = "black") +
annotate("text", x = 1.5, y = 11, label = "Wilcoxon test
p-value = .031 (V=21)")
# boxplot s8reln2
fl <- ggplot(mod_stim, aes(condition, s8reln2))
fl_p <- fl + geom_boxplot(outlier.shape = NA) +
geom_jitter(width = 0.1) +
labs(
title = "Final Lengthening on name2 by condition for model stimuli",
) +
theme_bw() +
theme(axis.text.x = element_text(size=14),
axis.text.y = element_text(size=14),
axis.title=element_text(size=14)) +
annotate("segment", x = 1, xend = 2, y = 49, yend = 49,
colour = "black") +
annotate("text", x = 1.5, y = 48, label = "Wilcoxon test
p-value = .031 (V=21)")
# boxplot pause3
p3 <- ggplot(mod_stim, aes(condition, pause3))
p3_p <- p3 + geom_boxplot(outlier.shape = NA) +
geom_jitter(width = 0.1)+
labs(
title = "Pause after name2 by condition for model stimuli",
) +
theme_bw() +
theme(axis.text.x = element_text(size=14),
axis.text.y = element_text(size=14),
axis.title=element_text(size=14)) +
annotate("segment", x = 1, xend = 2, y = 24, yend = 24,
colour = "black") +
annotate("text", x = 1.5, y = 23, label = "Wilcoxon test
p-value = .031 (V=21)")
## plot grids =================================
require(ggpubr)
mod_stim_cues <- ggarrange(f0_p,fl_p,p3_p, ncol=3, nrow=1)
mod_stim_cues
# ggsave("mod_stim_cue.jpg", plot = mod_stim_cues, width = 18, height = 10, dpi = 300)14 Session information
sessionInfo(package = NULL)## R version 4.2.1 (2022-06-23)
## Platform: aarch64-apple-darwin20 (64-bit)
## Running under: macOS Monterey 12.5
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.2-arm64/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.2-arm64/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] flextable_0.7.2 ggpubr_0.4.0 lmerTest_3.1-3 lme4_1.1-30
## [5] Matrix_1.4-1 forcats_0.5.1 stringr_1.4.0 dplyr_1.0.9
## [9] purrr_0.3.4 readr_2.1.2 tidyr_1.2.0 tibble_3.1.8
## [13] ggplot2_3.3.6 tidyverse_1.3.2 plyr_1.8.7
##
## loaded via a namespace (and not attached):
## [1] TH.data_1.1-1 googledrive_2.0.0 minqa_1.2.4
## [4] colorspace_2.0-3 ggsignif_0.6.3 ellipsis_0.3.2
## [7] estimability_1.4 base64enc_0.1-3 fs_1.5.2
## [10] rstudioapi_0.13 farver_2.1.1 fansi_1.0.3
## [13] mvtnorm_1.1-3 lubridate_1.8.0 xml2_1.3.3
## [16] codetools_0.2-18 splines_4.2.1 cachem_1.0.6
## [19] knitr_1.39 jsonlite_1.8.0 nloptr_2.0.3
## [22] pbkrtest_0.5.1 broom_1.0.0 dbplyr_2.2.1
## [25] compiler_4.2.1 httr_1.4.3 emmeans_1.7.5
## [28] backports_1.4.1 assertthat_0.2.1 fastmap_1.1.0
## [31] gargle_1.2.0 cli_3.3.0 htmltools_0.5.3
## [34] tools_4.2.1 coda_0.19-4 gtable_0.3.0
## [37] glue_1.6.2 Rcpp_1.0.9 carData_3.0-5
## [40] cellranger_1.1.0 jquerylib_0.1.4 vctrs_0.4.1
## [43] nlme_3.1-157 xfun_0.31 rvest_1.0.2
## [46] lifecycle_1.0.1 rstatix_0.7.0 googlesheets4_1.0.0
## [49] MASS_7.3-57 zoo_1.8-10 scales_1.2.0
## [52] hms_1.1.1 parallel_4.2.1 sandwich_3.0-2
## [55] RColorBrewer_1.1-3 yaml_2.3.5 gridExtra_2.3
## [58] gdtools_0.2.4 sass_0.4.2 stringi_1.7.8
## [61] highr_0.9 zip_2.2.0 boot_1.3-28
## [64] systemfonts_1.0.4 rlang_1.0.4 pkgconfig_2.0.3
## [67] evaluate_0.15 lattice_0.20-45 labeling_0.4.2
## [70] cowplot_1.1.1 tidyselect_1.1.2 magrittr_2.0.3
## [73] R6_2.5.1 generics_0.1.3 multcomp_1.4-19
## [76] DBI_1.1.3 pillar_1.8.0 haven_2.5.0
## [79] withr_2.5.0 mgcv_1.8-40 survival_3.3-1
## [82] abind_1.4-5 modelr_0.1.8 crayon_1.5.1
## [85] car_3.1-0 uuid_1.1-0 utf8_1.2.2
## [88] officer_0.4.3 tzdb_0.3.0 rmarkdown_2.14
## [91] grid_4.2.1 readxl_1.4.0 data.table_1.14.2
## [94] reprex_2.0.1 digest_0.6.29 xtable_1.8-4
## [97] numDeriv_2016.8-1.1 munsell_0.5.0 bslib_0.4.0rm(list=ls())